/[escript]/trunk/paso/src/Solvers/Solver_applyBlockDiagonalMatrix.c
ViewVC logotype

Diff of /trunk/paso/src/Solvers/Solver_applyBlockDiagonalMatrix.c

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 494 by jgs, Wed Nov 9 02:02:19 2005 UTC revision 495 by gross, Mon Feb 6 06:32:06 2006 UTC
# Line 21  Line 21 
21    
22  void Paso_Solver_applyBlockDiagonalMatrix(dim_t n_block,dim_t n,double* D,index_t* pivot,double* x,double* b) {  void Paso_Solver_applyBlockDiagonalMatrix(dim_t n_block,dim_t n,double* D,index_t* pivot,double* x,double* b) {
23       dim_t i;       dim_t i;
24         register dim_t i3,i9;
25         register double b0,b1,b2,D00,D10,D20,D01,D11,D21,D02,D12,D22;
26    
27       if (n_block==1) {       if (n_block==1) {
28           #pragma omp for private(i) schedule(static)           #pragma omp for private(i) schedule(static)
29           for (i=0;i<n;i++) {           for (i=0;i<n;++i) {
30              x[i]=D[i]*b[i];              x[i]=D[i]*b[i];
31           }           }
32       } else if (n_block==2) {       } else if (n_block==2) {
33           #pragma omp for private(i) schedule(static)           #pragma omp for private(i,b0,b1,D00,D10,D01,D11,i3,i9) schedule(static)
34           for (i=0;i<n;i++) {           for (i=0;i<n;++i) {
35              x[2*i  ]=D[4*i  ]*b[2*i]+D[4*i+2]*b[2*i+1];              i3=2*i;
36              x[2*i+1]=D[4*i+1]*b[2*i]+D[4*i+3]*b[2*i+1];              i9=4*i;
37                b0=b[i3];
38                b1=b[i3+1];
39                D00=D[i9  ];
40                D10=D[i9+1];
41                D01=D[i9+2];
42                D11=D[i9+3];
43                x[i3  ]=D00*b0+D01*b1;
44                x[i3+1]=D10*b0+D11*b1;
45           }           }
46       } else if (n_block==3) {       } else if (n_block==3) {
47           #pragma omp for private(i) schedule(static)           #pragma omp for private(i,b0,b1,b2,D00,D10,D20,D01,D11,D21,D02,D12,D22,i3,i9) schedule(static)
48           for (i=0;i<n;i++) {           for (i=0;i<n;++i) {
49              x[3*i  ]=D[9*i  ]*b[3*i]+D[9*i+3]*b[3*i+1]+D[9*i+6]*b[3*i+2];              i3=3*i;
50              x[3*i+1]=D[9*i+1]*b[3*i]+D[9*i+4]*b[3*i+1]+D[9*i+7]*b[3*i+2];              i9=9*i;
51              x[3*i+2]=D[9*i+2]*b[3*i]+D[9*i+5]*b[3*i+1]+D[9*i+8]*b[3*i+2];              b0=b[i3];
52                b1=b[i3+1];
53                b2=b[i3+2];
54                D00=D[i9  ];
55                D10=D[i9+1];
56                D20=D[i9+2];
57                D01=D[i9+3];
58                D11=D[i9+4];
59                D21=D[i9+5];
60                D02=D[i9+6];
61                D12=D[i9+7];
62                D22=D[i9+8];
63                x[i3  ]=D00*b0+D01*b1+D02*b2;
64                x[i3+1]=D10*b0+D11*b1+D12*b2;
65                x[i3+2]=D20*b0+D21*b1+D22*b2;
66           }           }
67       }       }
68       return;       return;

Legend:
Removed from v.494  
changed lines
  Added in v.495

  ViewVC Help
Powered by ViewVC 1.1.26