Loading...
Searching...
No Matches
blas_level_3.c
Go to the documentation of this file.
/*****************************************************************************
*
* MODULE: Grass numerical math interface
* AUTHOR(S): Soeren Gebbert, Berlin (GER) Dec 2006
* soerengebbert <at> googlemail <dot> com
*
* PURPOSE: grass blas implementation
* part of the gmath library
*
* COPYRIGHT: (C) 2010 by the GRASS Development Team
*
* This program is free software under the GNU General Public
* License (>=v2). Read the file COPYING that comes with GRASS
* for details.
*
*****************************************************************************/
#include <math.h>
#include <unistd.h>
#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#include <grass/gis.h>
#include <grass/gmath.h>
/*!
* \brief Add two matrices and scale matrix A with the scalar a
*
* \f[ {\bf C} = a {\bf A} + {\bf B} \f]
*
* In case B == NULL, matrix A will be scaled by scalar a. \n
* In case a == 1.0, a simple matrix addition is performed. \n
* In case a == -1.0 matrix A is subtracted from matrix B. \n
* The result is written into matrix C.
*
*
* This function is multi-threaded with OpenMP and can be called within a
* parallel OpenMP region.
*
* \param A (double **)
* \param B (double **) if NULL, matrix A is scaled by scalar a only
* \param a (double)
* \param C (double **)
* \param rows (int)
* \param cols (int)
* \return (void)
*
* */
void G_math_d_aA_B(double **A, double **B, double a, double **C, int rows,
int cols)
{
int i, j;
/*If B is null, scale the matrix A with th scalar a */
if (B == NULL) {
#pragma omp for schedule(static) private(i, j)
for (i = rows - 1; i >= 0; i--)
for (j = cols - 1; j >= 0; j--)
C[i][j] = a * A[i][j];
return;
}
/*select special cases */
if (a == 1.0) {
#pragma omp for schedule(static) private(i, j)
for (i = rows - 1; i >= 0; i--)
for (j = cols - 1; j >= 0; j--)
C[i][j] = A[i][j] + B[i][j];
}
else if (a == -1.0) {
#pragma omp for schedule(static) private(i, j)
for (i = rows - 1; i >= 0; i--)
for (j = cols - 1; j >= 0; j--)
C[i][j] = B[i][j] - A[i][j];
}
else {
#pragma omp for schedule(static) private(i, j)
for (i = rows - 1; i >= 0; i--)
for (j = cols - 1; j >= 0; j--)
C[i][j] = a * A[i][j] + B[i][j];
}
return;
}
/*!
* \brief Add two matrices and scale matrix A with the scalar a
*
* \f[ {\bf C} = a {\bf A} + {\bf B} \f]
*
* In case B == NULL, matrix A will be scaled by scalar a. \n
* In case a == 1.0, a simple matrix addition is performed. \n
* In case a == -1.0 matrix A is subtracted from matrix B. \n
* The result is written into matrix C.
*
*
*
* This function is multi-threaded with OpenMP and can be called within a
parallel OpenMP region.
*
* \param A (float **)
* \param B (float **) if NULL, matrix A is scaled by scalar a only
* \param a (float)
* \param C (float **)
* \param rows (int)
* \param cols (int)
* \return (void)
*
* */
void G_math_f_aA_B(float **A, float **B, float a, float **C, int rows, int cols)
{
int i, j;
/*If B is null, scale the matrix A with th scalar a */
if (B == NULL) {
#pragma omp for schedule(static) private(i, j)
for (i = rows - 1; i >= 0; i--)
for (j = cols - 1; j >= 0; j--)
C[i][j] = a * A[i][j];
return;
}
/*select special cases */
if (a == 1.0) {
#pragma omp for schedule(static) private(i, j)
for (i = rows - 1; i >= 0; i--)
for (j = cols - 1; j >= 0; j--)
C[i][j] = A[i][j] + B[i][j];
}
else if (a == -1.0) {
#pragma omp for schedule(static) private(i, j)
for (i = rows - 1; i >= 0; i--)
for (j = cols - 1; j >= 0; j--)
C[i][j] = B[i][j] - A[i][j];
}
else {
#pragma omp for schedule(static) private(i, j)
for (i = rows - 1; i >= 0; i--)
for (j = cols - 1; j >= 0; j--)
C[i][j] = a * A[i][j] + B[i][j];
}
return;
}
/*!
* \brief Matrix multiplication
*
* \f[ {\bf C} = {\bf A}{\bf B} \f]
*
* The result is written into matrix C.
*
* A must be of size rows_A * cols_A
* B must be of size rows_B * cols_B with rows_B == cols_A
* C must be of size rows_A * cols_B
*
*
* This function is multi-threaded with OpenMP and can be called within a
* parallel OpenMP region.
*
* \param A (double **)
* \param B (double **)
* \param C (double **)
* \param rows_A (int)
* \param cols_A (int)
* \param cols_B (int)
* \return (void)
*
* */
void G_math_d_AB(double **A, double **B, double **C, int rows_A, int cols_A,
int cols_B)
{
int i, j, k;
#pragma omp for schedule(static) private(i, j, k)
for (i = 0; i < rows_A; i++) {
for (j = 0; j < cols_B; j++) {
C[i][j] = 0.0;
for (k = cols_A - 1; k >= 0; k--) {
C[i][j] += A[i][k] * B[k][j];
}
}
}
return;
}
/*!
* \brief Matrix multiplication
*
* \f[ {\bf C} = {\bf A}{\bf B} \f]
*
* The result is written into matrix C.
*
* A must be of size rows_A * cols_A
* B must be of size rows_B * cols_B with rows_B == cols_A
* C must be of size rows_A * cols_B
*
*
* This function is multi-threaded with OpenMP and can be called within a
* parallel OpenMP region.
*
* \param A (float **)
* \param B (float **)
* \param C (float **)
* \param rows_A (int)
* \param cols_A (int)
* \param cols_B (int)
* \return (void)
*
* */
void G_math_f_AB(float **A, float **B, float **C, int rows_A, int cols_A,
int cols_B)
{
int i, j, k;
#pragma omp for schedule(static) private(i, j, k)
for (i = 0; i < rows_A; i++) {
for (j = 0; j < cols_B; j++) {
C[i][j] = 0.0;
for (k = cols_A - 1; k >= 0; k--) {
C[i][j] += A[i][k] * B[k][j];
}
}
}
return;
}
G_math_f_aA_B
void G_math_f_aA_B(float **A, float **B, float a, float **C, int rows, int cols)
Add two matrices and scale matrix A with the scalar a.
Definition blas_level_3.c:112
G_math_f_AB
void G_math_f_AB(float **A, float **B, float **C, int rows_A, int cols_A, int cols_B)
Matrix multiplication.
Definition blas_level_3.c:214
G_math_d_aA_B
void G_math_d_aA_B(double **A, double **B, double a, double **C, int rows, int cols)
Add two matrices and scale matrix A with the scalar a.
Definition blas_level_3.c:49
G_math_d_AB
void G_math_d_AB(double **A, double **B, double **C, int rows_A, int cols_A, int cols_B)
Matrix multiplication.
Definition blas_level_3.c:172
NULL
#define NULL
Definition ccmath.h:32