[go: up one dir, main page]
Menu ▾ ▴

[r7]: / src / assembLoop.c  Maximize  Restore  History

Download this file

53 lines (43 with data), 1.4 kB 

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
#include "mex.h"

#include <math.h>



/* The computational routine */

void assembLoop(double* bDet, double *tri, double *f, double *b, mwSize n)

{

    double S4[3][3] = { {2, 1, 1},

                        {1, 2, 1},

                        {1, 1, 2}};



    mwSize i,j,m;

    double temp, myDet; 

                        

    for (i = 0; i < n; i++) {

      myDet = fabs(bDet[i]);

      for (m = 0; m < 3; m++) {

        temp = 0;

        for (j = 0; j < 3; j++) {

          temp += f[(int)tri[i + n*j]-1] * S4[m][j] / 24;

        }

        b[(int)tri[i + n*m]-1] += myDet * temp;

      }

    }

}



void mexFunction( int nlhs, mxArray *plhs[],

                  int nrhs, const mxArray *prhs[])

{

    double* bDet;

    double* tri;

    double* f;

    size_t columnsF,rowsF,columnsT,rowsT;

    double* b;



    /* Inputs parsen und analysisieren */

    bDet = mxGetPr(prhs[0]);

    tri  =  mxGetPr(prhs[1]);

    f    = mxGetPr(prhs[2]);



    rowsF    = (int) mxGetM(prhs[2]);

    columnsF = (int) mxGetN(prhs[2]);



    rowsT    = (int) mxGetM(prhs[1]);

    columnsT = (int) mxGetN(prhs[1]);



    /* Output Erzeugen */

    plhs[0] = mxCreateDoubleMatrix((mwSize)rowsF, (mwSize)columnsF, mxREAL);

    b       = mxGetPr(plhs[0]);



    /* call the computational routine */

    assembLoop(bDet,tri,f,b,(mwSize)rowsT);

}