Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
72 changes: 70 additions & 2 deletions pygrt/C_extension/include/grt/common/matrix.h
Original file line number Diff line number Diff line change
Expand Up @@ -3,7 +3,7 @@
* @author Zhu Dengda (zhudengda@mail.iggcas.ac.cn)
* @date 2024-07-24
*
* 2x2小矩阵的加、减、乘、除、求逆等操作,由于均为小型数组操作,所有函数均为内联函数
* 小矩阵的加、减、乘、除、求逆等操作,由于均为小型数组操作,所有函数均为内联函数
*/

#pragma once
Expand Down Expand Up @@ -142,7 +142,7 @@ inline GCC_ALWAYS_INLINE void grt_cmat2x2_assign(const MYCOMPLEX M1[2][2], MYCOM
}

/**
* 计算nxn复矩阵的积(小矩阵)(最暴力的方式)
* 计算mxn复矩阵的积(小矩阵)(最暴力的方式)
*
* @param[in] m1 M1矩阵行数
* @param[in] n1 M1矩阵列数
Expand Down Expand Up @@ -171,6 +171,31 @@ inline GCC_ALWAYS_INLINE void grt_cmatmxn_mul(MYINT m1, MYINT n1, MYINT p1, cons
}
}

/**
* 计算mxn复矩阵的转置矩阵(不共轭)
*
* @param[in] m1 M1矩阵行数
* @param[in] n1 M1矩阵列数
* @param[in] M1 M1矩阵
* @param[out] M2 M2矩阵 (M1^T)
*/
inline GCC_ALWAYS_INLINE void grt_cmatmxn_transpose(MYINT m1, MYINT n1, const MYCOMPLEX M1[m1][n1], MYCOMPLEX M2[n1][m1]){
MYINT m, n;
MYCOMPLEX M0[n1][m1];
for(m=0; m<m1; ++m){
for(n=0; n<n1; ++n){
M0[n][m] = M1[m][n];
}
}

// memcpy(M, M0, sizeof(MYCOMPLEX)*m1*p1);
for(m=0; m<m1; ++m){
for(n=0; n<n1; ++n){
M2[n][m] = M0[n][m];
}
}
}

/**
* 从M1大矩阵中划分Q子矩阵
*
Expand All @@ -191,6 +216,28 @@ inline GCC_ALWAYS_INLINE void grt_cmatmxn_block(MYINT m1, MYINT n1, const MYCOMP
}
}


/**
* 将小矩阵Q填充到M1大矩阵中
*
* @param[in] m1 M1矩阵行数
* @param[in] n1 M1矩阵列数
* @param[in] M1 M1矩阵
* @param[in] im 子矩阵起始行索引
* @param[in] in 子矩阵起始列索引
* @param[in] lm 子矩阵行数
* @param[in] ln 子矩阵列数
* @param[out] Q 子矩阵
*/
inline GCC_ALWAYS_INLINE void grt_cmatmxn_block_assign(MYINT m1, MYINT n1, MYCOMPLEX M[m1][n1], MYINT im, MYINT in, MYINT lm, MYINT ln, const MYCOMPLEX Q[lm][ln]){
for(MYINT m=0; m<lm; ++m){
for(MYINT n=0; n<ln; ++n){
M[im+m][in+n] = Q[m][n];
}
}
}


/**
* 打印矩阵
*
Expand All @@ -206,4 +253,25 @@ inline GCC_ALWAYS_INLINE void grt_cmatmxn_print(MYINT m1, MYINT n1, const MYCOMP
}
fprintf(stderr, "\n");
}
}


/**
* 计算mxn复矩阵和nx1的复向量的积
*
* @param[in] M1 矩阵1
* @param[in] M2 向量
* @param[out] M 积矩阵, M1 * M2
*/
inline GCC_ALWAYS_INLINE void grt_cmatmx1_mul(MYINT m1, MYINT n1, const MYCOMPLEX M1[m1][n1], const MYCOMPLEX M2[n1], MYCOMPLEX M[n1]){
MYCOMPLEX M0[n1];
for(MYINT i=0; i<m1; ++i){
M0[i] = 0.0;
for(MYINT j=0; j<n1; ++j){
M0[i] += M1[i][j] * M2[j];
}
}
for(MYINT i=0; i<n1; ++i){
M[i] = M0[i];
}
}
87 changes: 74 additions & 13 deletions pygrt/C_extension/include/grt/common/search.h
Original file line number Diff line number Diff line change
Expand Up @@ -11,11 +11,16 @@
#include <stdbool.h>
#include "grt/common/const.h"

// 定义 X 宏,为多个类型定义查找函数
#define __FOR_EACH_TYPE \
X(MYINT) X(MYREAL) X(float) X(double)


/**
* 该函数对输入的整数数组进行线性搜索,找到目标值时返回其索引。
* 该函数对输入数组进行线性搜索,找到目标值时返回其索引。
* 如果目标值在数组中未找到,则返回 -1。
*
* @param[in] array 要搜索的整数数组
* @param[in] array 要搜索的数组
* @param[in] size 数组的大小(元素个数)。
* @param[in] target 要查找的目标值。
*
Expand All @@ -24,13 +29,19 @@
* @note 如果数组中存在多个目标值,该函数返回第一个匹配的索引。
*
*/
MYINT grt_findElement_MYINT(const MYINT array[], MYINT size, MYINT target);
#define X(T) \
MYINT grt_findElement_##T(const T *array, MYINT size, T target);

__FOR_EACH_TYPE
#undef X



/**
* 搜索浮点数数组中最接近目标值且小于目标值的索引
* 搜索数组中最接近目标值且小于目标值的索引
* 如果目标值在数组中未找到,则返回 -1。
*
* @param[in] array 要搜索的浮点数数组
* @param[in] array 要搜索的数组
* @param[in] size 数组的大小(元素个数)。
* @param[in] target 要查找的目标值。
*
Expand All @@ -39,12 +50,18 @@ MYINT grt_findElement_MYINT(const MYINT array[], MYINT size, MYINT target);
* @note 如果数组中存在多个目标值,该函数返回第一个匹配的索引。
*
*/
MYINT grt_findLessEqualClosest_MYREAL(const MYREAL array[], MYINT size, MYREAL target);
#define X(T) \
MYINT grt_findLessEqualClosest_##T(const T *array, MYINT size, T target);

__FOR_EACH_TYPE
#undef X



/**
* 搜索浮点数数组中最接近目标值的索引
* 搜索数组中最接近目标值的索引
*
* @param[in] array 要搜索的浮点数数组
* @param[in] array 要搜索的数组
* @param[in] size 数组的大小(元素个数)。
* @param[in] target 要查找的目标值。
*
Expand All @@ -53,18 +70,62 @@ MYINT grt_findLessEqualClosest_MYREAL(const MYREAL array[], MYINT size, MYREAL t
* @note 如果数组中存在多个目标值,该函数返回第一个匹配的索引。
*
*/
MYINT grt_findClosest_MYREAL(const MYREAL array[], MYINT size, MYREAL target);
#define X(T) \
MYINT grt_findClosest_##T(const T *array, MYINT size, T target);

__FOR_EACH_TYPE
#undef X


/**
* 搜索浮点数数组的最大或最小值,返回其索引。
* 搜索数组的最值,返回其索引。
*
* @param[in] array 要搜索的浮点数数组
* @param[in] array 要搜索的数组
* @param[in] size 数组的大小(元素个数)。
* @param[in] isMax 是否要找最大值,否则找最小值。
*
* @return idx 目标值的索引。
*
* @note 如果数组中存在相同最值,该函数返回第一个匹配的索引。
*
*/
MYINT grt_findMinMax_MYREAL(const MYREAL array[], MYINT size, bool isMax);
#define X(T) \
MYINT grt_findMin_##T(const T *array, MYINT size);\
MYINT grt_findMax_##T(const T *array, MYINT size);\

__FOR_EACH_TYPE
#undef X


/**
* 比较函数
*
* @param[in] a 元素 a 地址
* @param[in] b 元素 b 地址
*
* @return flag 比较结果,(1) a > b, (0) a == b, (-1) a < b
*/
#define X(T) \
int grt_compare_##T(const void *a, const void *b);

__FOR_EACH_TYPE
#undef X
#undef __FOR_EACH_TYPE



/**
* 在有序数组中插入元素,元素类型和数组类型需匹配
*
* @param[in,out] arr 有序数组地址
* @param[in,out] size 数组大小地址
* @param[in] capacity 数组最大容量
* @param[in] target 元素地址
* @param[in] elementSize 元素和数组内元素的字节长度
* @param[in] ascending 升序(true) 或 降序(false)
* @param[in] compare 比较函数
*
* @return pos 插入位置的索引
*/
MYINT grt_insertOrdered(
void *arr, MYINT *size, MYINT capacity, const void *target, size_t elementSize, bool ascending,
int (*compare)(const void *, const void *));
Loading