aboutsummaryrefslogtreecommitdiff
path: root/libcrystfel/src
diff options
context:
space:
mode:
Diffstat (limited to 'libcrystfel/src')
-rw-r--r--libcrystfel/src/cell-utils.c46
-rw-r--r--libcrystfel/src/integer_matrix.c4
-rw-r--r--libcrystfel/src/integer_matrix.h4
-rw-r--r--libcrystfel/src/rational.c61
-rw-r--r--libcrystfel/src/rational.h11
-rw-r--r--libcrystfel/src/symmetry.c10
6 files changed, 72 insertions, 64 deletions
diff --git a/libcrystfel/src/cell-utils.c b/libcrystfel/src/cell-utils.c
index dbc845bc..de8fe7ae 100644
--- a/libcrystfel/src/cell-utils.c
+++ b/libcrystfel/src/cell-utils.c
@@ -1612,7 +1612,6 @@ int compare_derivative_cell_parameters(UnitCell *cell_in, UnitCell *reference_in
Rational *cand_c;
int ncand_a, ncand_b, ncand_c;
int ia, ib;
- RationalMatrix *MCB;
RationalMatrix *CiAMCB;
double min_dist = +INFINITY;
@@ -1647,8 +1646,6 @@ int compare_derivative_cell_parameters(UnitCell *cell_in, UnitCell *reference_in
}
M = rtnl_mtx_new(3, 3);
- MCB = rtnl_mtx_new(3, 3);
- CiAMCB = rtnl_mtx_new(3, 3);
for ( ia=0; ia<ncand_a; ia++ ) {
for ( ib=0; ib<ncand_b; ib++ ) {
@@ -1656,6 +1653,7 @@ int compare_derivative_cell_parameters(UnitCell *cell_in, UnitCell *reference_in
double at, bt, ct, alt, bet, gat;
double dist;
int ic = 0;
+ RationalMatrix *MCB;
/* Form the matrix using the first candidate for c */
rtnl_mtx_set(M, 0, 0, cand_a[3*ia+0]);
@@ -1710,8 +1708,9 @@ int compare_derivative_cell_parameters(UnitCell *cell_in, UnitCell *reference_in
dist = g6_distance(test, reference);
if ( dist < min_dist ) {
min_dist = dist;
- rtnl_mtx_mtxmult(M, CB, MCB);
- rtnl_mtx_mtxmult(CiA, MCB, CiAMCB);
+ MCB = rtnlmtx_times_rtnlmtx(M, CB);
+ CiAMCB = rtnlmtx_times_rtnlmtx(CiA, MCB);
+ rtnl_mtx_free(MCB);
}
cell_free(test);
@@ -1721,7 +1720,6 @@ int compare_derivative_cell_parameters(UnitCell *cell_in, UnitCell *reference_in
}
rtnl_mtx_free(M);
- rtnl_mtx_free(MCB);
free(cand_a);
free(cand_b);
free(cand_c);
@@ -1817,7 +1815,7 @@ static void debug_lattice(struct g6 g, double eps, int step)
static void mult_in_place(IntegerMatrix *T, IntegerMatrix *M)
{
int i, j;
- IntegerMatrix *tmp = intmat_intmat_mult(T, M);
+ IntegerMatrix *tmp = intmat_times_intmat(T, M);
for ( i=0; i<3; i++ ) {
for ( j=0; j<3; j++ ) {
intmat_set(T, i, j, intmat_get(tmp, i, j));
@@ -1827,40 +1825,6 @@ static void mult_in_place(IntegerMatrix *T, IntegerMatrix *M)
}
-/* Replace the elements of T with those of T*M */
-static void rtnl_mult_in_place(RationalMatrix *T, RationalMatrix *M)
-{
- int i, j;
- RationalMatrix *tmp;
-
- tmp = rtnl_mtx_new(3, 3);
- rtnl_mtx_mtxmult(T, M, tmp);
- for ( i=0; i<3; i++ ) {
- for ( j=0; j<3; j++ ) {
- rtnl_mtx_set(T, i, j, rtnl_mtx_get(tmp, i, j));
- }
- }
- rtnl_mtx_free(tmp);
-}
-
-
-/* Replace the elements of T with those of T*M */
-static void rtnl_int_mult_in_place(RationalMatrix *T, IntegerMatrix *M)
-{
- int i, j;
- RationalMatrix *tmp;
-
- tmp = rtnl_mtx_new(3, 3);
- rtnl_mtx_intmatmult(T, M, tmp);
- for ( i=0; i<3; i++ ) {
- for ( j=0; j<3; j++ ) {
- rtnl_mtx_set(T, i, j, rtnl_mtx_get(tmp, i, j));
- }
- }
- rtnl_mtx_free(tmp);
-}
-
-
/* Cell volume from G6 components */
static double g6_volume(struct g6 g)
{
diff --git a/libcrystfel/src/integer_matrix.c b/libcrystfel/src/integer_matrix.c
index c633a620..16aff5bc 100644
--- a/libcrystfel/src/integer_matrix.c
+++ b/libcrystfel/src/integer_matrix.c
@@ -218,8 +218,8 @@ signed int *transform_indices(const IntegerMatrix *P, const signed int *hkl)
* \returns A newly allocated \ref IntegerMatrix containing the answer, or NULL on
* error.
**/
-IntegerMatrix *intmat_intmat_mult(const IntegerMatrix *a,
- const IntegerMatrix *b)
+IntegerMatrix *intmat_times_intmat(const IntegerMatrix *a,
+ const IntegerMatrix *b)
{
unsigned int i, j;
IntegerMatrix *ans;
diff --git a/libcrystfel/src/integer_matrix.h b/libcrystfel/src/integer_matrix.h
index 5a365f0a..25c0fb99 100644
--- a/libcrystfel/src/integer_matrix.h
+++ b/libcrystfel/src/integer_matrix.h
@@ -74,8 +74,8 @@ extern IntegerMatrix *intmat_create_3x3(signed int m11, signed int m12, signed i
extern signed int *transform_indices(const IntegerMatrix *P, const signed int *hkl);
/* Matrix-matrix multiplication */
-extern IntegerMatrix *intmat_intmat_mult(const IntegerMatrix *a,
- const IntegerMatrix *b);
+extern IntegerMatrix *intmat_times_intmat(const IntegerMatrix *a,
+ const IntegerMatrix *b);
/* Inverse */
extern IntegerMatrix *intmat_inverse(const IntegerMatrix *m);
diff --git a/libcrystfel/src/rational.c b/libcrystfel/src/rational.c
index f8d3cb68..d59da59c 100644
--- a/libcrystfel/src/rational.c
+++ b/libcrystfel/src/rational.c
@@ -528,14 +528,18 @@ void rtnl_mtx_print(const RationalMatrix *m)
}
-void rtnl_mtx_mtxmult(const RationalMatrix *A, const RationalMatrix *B,
- RationalMatrix *ans)
+RationalMatrix *rtnlmtx_times_intmat(const RationalMatrix *A,
+ const IntegerMatrix *B)
{
int i, j;
+ RationalMatrix *ans;
+ unsigned int B_rows, B_cols;
- assert(ans->cols == B->cols);
- assert(ans->rows == A->rows);
- assert(A->cols == B->rows);
+ intmat_size(B, &B_rows, &B_cols);
+ assert(A->cols == B_rows);
+
+ ans = rtnl_mtx_new(A->rows, B_cols);
+ if ( ans == NULL ) return NULL;
for ( i=0; i<ans->rows; i++ ) {
for ( j=0; j<ans->cols; j++ ) {
@@ -544,21 +548,27 @@ void rtnl_mtx_mtxmult(const RationalMatrix *A, const RationalMatrix *B,
for ( k=0; k<A->rows; k++ ) {
Rational add;
add = rtnl_mul(rtnl_mtx_get(A, i, k),
- rtnl_mtx_get(B, k, j));
+ rtnl(intmat_get(B, k, j), 1));
sum = rtnl_add(sum, add);
}
rtnl_mtx_set(ans, i, j, sum);
}
}
+
+ return ans;
}
-void rtnl_mtx_intmatmult(const RationalMatrix *A, const IntegerMatrix *B,
- RationalMatrix *ans)
+RationalMatrix *rtnlmtx_times_rtnlmtx(const RationalMatrix *A,
+ const RationalMatrix *B)
{
int i, j;
+ RationalMatrix *ans;
+
+ assert(A->cols == B->rows);
- assert(ans->rows == A->rows);
+ ans = rtnl_mtx_new(A->rows, B->cols);
+ if ( ans == NULL ) return NULL;
for ( i=0; i<ans->rows; i++ ) {
for ( j=0; j<ans->cols; j++ ) {
@@ -567,12 +577,43 @@ void rtnl_mtx_intmatmult(const RationalMatrix *A, const IntegerMatrix *B,
for ( k=0; k<A->rows; k++ ) {
Rational add;
add = rtnl_mul(rtnl_mtx_get(A, i, k),
- rtnl(intmat_get(B, k, j), 1));
+ rtnl_mtx_get(B, k, j));
+ sum = rtnl_add(sum, add);
+ }
+ rtnl_mtx_set(ans, i, j, sum);
+ }
+ }
+
+ return ans;
+}
+
+
+RationalMatrix *intmat_times_rtnlmtx(const IntegerMatrix *A, const RationalMatrix *B)
+{
+ int i, j;
+ RationalMatrix *ans;
+ unsigned int A_rows, A_cols;
+
+ intmat_size(A, &A_rows, &A_cols);
+
+ ans = rtnl_mtx_new(A_rows, B->cols);
+ if ( ans == NULL ) return NULL;
+
+ for ( i=0; i<ans->rows; i++ ) {
+ for ( j=0; j<ans->cols; j++ ) {
+ int k;
+ Rational sum = rtnl_zero();
+ for ( k=0; k<A_rows; k++ ) {
+ Rational add;
+ add = rtnl_mul(rtnl(intmat_get(A, i, k), 1),
+ rtnl_mtx_get(B, k, j));
sum = rtnl_add(sum, add);
}
rtnl_mtx_set(ans, i, j, sum);
}
}
+
+ return ans;
}
diff --git a/libcrystfel/src/rational.h b/libcrystfel/src/rational.h
index 012b929e..880bdbeb 100644
--- a/libcrystfel/src/rational.h
+++ b/libcrystfel/src/rational.h
@@ -90,10 +90,13 @@ extern RationalMatrix *rtnl_mtx_from_intmat(const IntegerMatrix *m);
extern RationalMatrix *rtnl_mtx_identity(int rows);
extern IntegerMatrix *intmat_from_rtnl_mtx(const RationalMatrix *m);
extern void rtnl_mtx_free(RationalMatrix *mtx);
-extern void rtnl_mtx_mtxmult(const RationalMatrix *A, const RationalMatrix *B,
- RationalMatrix *ans);
-extern void rtnl_mtx_intmatmult(const RationalMatrix *A, const IntegerMatrix *B,
- RationalMatrix *ans);
+extern RationalMatrix *rtnlmtx_times_rtnlmtx(const RationalMatrix *A,
+ const RationalMatrix *B);
+extern RationalMatrix *rtnlmtx_times_intmat(const RationalMatrix *A,
+ const IntegerMatrix *B);
+extern RationalMatrix *intmat_times_rtnlmtx(const IntegerMatrix *a,
+ const RationalMatrix *b);
+
extern int transform_fractional_coords_rtnl(const RationalMatrix *P,
const Rational *ivec,
Rational *ans);
diff --git a/libcrystfel/src/symmetry.c b/libcrystfel/src/symmetry.c
index 6739daaf..e834e357 100644
--- a/libcrystfel/src/symmetry.c
+++ b/libcrystfel/src/symmetry.c
@@ -326,7 +326,7 @@ static void expand_ops(SymOpList *s)
int k, nk;
int found;
- m = intmat_intmat_mult(opi, opj);
+ m = intmat_times_intmat(opi, opj);
assert(m != NULL);
nk = num_ops(s);
@@ -380,13 +380,13 @@ static void transform_ops(SymOpList *s, IntegerMatrix *P)
IntegerMatrix *r, *f;
- r = intmat_intmat_mult(P, s->ops[i]);
+ r = intmat_times_intmat(P, s->ops[i]);
if ( r == NULL ) {
ERROR("Matrix multiplication failed.\n");
return;
}
- f = intmat_intmat_mult(r, Pi);
+ f = intmat_times_intmat(r, Pi);
if ( f == NULL ) {
ERROR("Matrix multiplication failed.\n");
return;
@@ -1340,7 +1340,7 @@ static int check_mult(const IntegerMatrix *ans,
int val;
IntegerMatrix *m;
- m = intmat_intmat_mult(a, b);
+ m = intmat_times_intmat(a, b);
assert(m != NULL);
val = intmat_equals(ans, m);
@@ -1403,7 +1403,7 @@ static int order(const IntegerMatrix *m)
IntegerMatrix *anew;
- anew = intmat_intmat_mult(m, a);
+ anew = intmat_times_intmat(m, a);
assert(anew != NULL);
intmat_free(a);
a = anew;
generated by cgit v1.2.3 (git 2.25.1) at 2024-11-05 15:34:07 +0000