/[escript]/trunk/finley/src/Assemble_LumpedSystem.c
ViewVC logotype

Diff of /trunk/finley/src/Assemble_LumpedSystem.c

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

trunk/esys2/finley/src/finleyC/Assemble_PDE.c revision 149 by jgs, Thu Sep 1 03:31:39 2005 UTC trunk/finley/src/finleyC/Assemble_PDE.c revision 155 by jgs, Wed Nov 9 02:02:19 2005 UTC
# Line 1  Line 1 
1  /* $Id$ */  /*
2     ******************************************************************************
3     *                                                                            *
4     *       COPYRIGHT  ACcESS 2003,2004,2005 -  All Rights Reserved              *
5     *                                                                            *
6     * This software is the property of ACcESS. No part of this code              *
7     * may be copied in any form or by any means without the expressed written    *
8     * consent of ACcESS.  Copying, use or modification of this software          *
9     * by any unauthorised person is illegal unless that person has a software    *
10     * license agreement with ACcESS.                                             *
11     *                                                                            *
12     ******************************************************************************
13    */
14    
15    
16  /**************************************************************/  /**************************************************************/
17    
# Line 28  Line 41 
41    
42  /**************************************************************/  /**************************************************************/
43    
44  /*   Copyrights by ACcESS Australia, 2003,2004 */  /*  Author: gross@access.edu.au */
45  /*   author: gross@access.edu.au */  /*  Version: $Id$ */
 /*   Version: $Id$ */  
46    
47  /**************************************************************/  /**************************************************************/
48    
 #include "escript/Data/DataC.h"  
 #include "Finley.h"  
49  #include "Assemble.h"  #include "Assemble.h"
 #include "NodeFile.h"  
 #include "ElementFile.h"  
50  #include "Util.h"  #include "Util.h"
51  #ifdef _OPENMP  #ifdef _OPENMP
52  #include <omp.h>  #include <omp.h>
# Line 47  Line 55 
55    
56  /**************************************************************/  /**************************************************************/
57    
58  void Finley_Assemble_PDE(Finley_NodeFile* nodes,Finley_ElementFile* elements,Finley_SystemMatrix* S, escriptDataC* F,  void Finley_Assemble_PDE(Finley_NodeFile* nodes,Finley_ElementFile* elements,Paso_SystemMatrix* S, escriptDataC* F,
59               escriptDataC* A, escriptDataC* B, escriptDataC* C, escriptDataC* D, escriptDataC* X, escriptDataC* Y ) {               escriptDataC* A, escriptDataC* B, escriptDataC* C, escriptDataC* D, escriptDataC* X, escriptDataC* Y ) {
60    
61      char error_msg[LenErrorMsg_MAX];
62    double *EM_S=NULL,*EM_F=NULL,*V=NULL,*dVdv=NULL,*dSdV=NULL,*Vol=NULL,*dvdV=NULL;    double *EM_S=NULL,*EM_F=NULL,*V=NULL,*dVdv=NULL,*dSdV=NULL,*Vol=NULL,*dvdV=NULL;
63    double time0;    double time0;
64    dim_t dimensions[ESCRIPT_MAX_DATA_RANK],e,q;    dim_t dimensions[ESCRIPT_MAX_DATA_RANK],e,q;
65    Assemble_Parameters p;    Assemble_Parameters p;
66    index_t *index_row=NULL,*index_col=NULL,color;    index_t *index_row=NULL,*index_col=NULL,color;
67      Finley_resetError();
68    
69    if (nodes==NULL || elements==NULL) return;    if (nodes==NULL || elements==NULL) return;
70    if (S==NULL && isEmpty(F)) return;    if (S==NULL && isEmpty(F)) return;
71    
72    /* set all parameters in p*/    /* set all parameters in p*/
73    Assemble_getAssembleParameters(nodes,elements,S,F,&p);    Assemble_getAssembleParameters(nodes,elements,S,F,&p);
74    if (Finley_ErrorCode!=NO_ERROR) return;    if (! Finley_noError()) return;
75    
76    /*  this function assumes NS=NN */    /*  this function assumes NS=NN */
77    if (p.NN!=p.NS) {    if (p.NN!=p.NS) {
78      Finley_ErrorCode=SYSTEM_ERROR;      Finley_setError(SYSTEM_ERROR,"__FILE__: for Finley_Assemble_PDE numNodes and numShapes have to be identical.");
     sprintf(Finley_ErrorMsg,"for Finley_Assemble_PDE numNodes and numShapes have to be identical.");  
79      return;      return;
80    }    }
81    if (p.numDim!=p.numElementDim) {    if (p.numDim!=p.numElementDim) {
82      Finley_ErrorCode=SYSTEM_ERROR;      Finley_setError(SYSTEM_ERROR,"__FILE__: Finley_Assemble_PDE accepts volume elements only.");
     sprintf(Finley_ErrorMsg,"Finley_Assemble_PDE accepts volume elements only.");  
83      return;      return;
84    }    }
85    /*  get a functionspace */    /*  get a functionspace */
# Line 87  void Finley_Assemble_PDE(Finley_NodeFile Line 95  void Finley_Assemble_PDE(Finley_NodeFile
95    /* check if all function spaces are the same */    /* check if all function spaces are the same */
96    
97    if (! functionSpaceTypeEqual(funcspace,A) ) {    if (! functionSpaceTypeEqual(funcspace,A) ) {
98          Finley_ErrorCode=TYPE_ERROR;          Finley_setError(TYPE_ERROR,"__FILE__: unexpected function space type for coefficient A");
         sprintf(Finley_ErrorMsg,"unexpected function space type for coefficient A");  
99    }    }
100    if (! functionSpaceTypeEqual(funcspace,B) ) {    if (! functionSpaceTypeEqual(funcspace,B) ) {
101          Finley_ErrorCode=TYPE_ERROR;          Finley_setError(TYPE_ERROR,"__FILE__: unexpected function space type for coefficient B");
         sprintf(Finley_ErrorMsg,"unexpected function space type for coefficient B");  
102    }    }
103    if (! functionSpaceTypeEqual(funcspace,C) ) {    if (! functionSpaceTypeEqual(funcspace,C) ) {
104          Finley_ErrorCode=TYPE_ERROR;          Finley_setError(TYPE_ERROR,"__FILE__: unexpected function space type for coefficient C");
         sprintf(Finley_ErrorMsg,"unexpected function space type for coefficient C");  
105    }    }
106    if (! functionSpaceTypeEqual(funcspace,D) ) {    if (! functionSpaceTypeEqual(funcspace,D) ) {
107          Finley_ErrorCode=TYPE_ERROR;          Finley_setError(TYPE_ERROR,"__FILE__: unexpected function space type for coefficient D");
         sprintf(Finley_ErrorMsg,"unexpected function space type for coefficient D");  
108    }    }
109    if (! functionSpaceTypeEqual(funcspace,X) ) {    if (! functionSpaceTypeEqual(funcspace,X) ) {
110          Finley_ErrorCode=TYPE_ERROR;          Finley_setError(TYPE_ERROR,"__FILE__: unexpected function space type for coefficient X");
         sprintf(Finley_ErrorMsg,"unexpected function space type for coefficient X");  
111    }    }
112    if (! functionSpaceTypeEqual(funcspace,Y) ) {    if (! functionSpaceTypeEqual(funcspace,Y) ) {
113          Finley_ErrorCode=TYPE_ERROR;          Finley_setError(TYPE_ERROR,"__FILE__: unexpected function space type for coefficient Y");
         sprintf(Finley_ErrorMsg,"unexpected function space type for coefficient Y");  
114    }    }
115    
116    /* check if all function spaces are the same */    /* check if all function spaces are the same */
117    
118    if (! numSamplesEqual(A,p.numQuad,elements->numElements) ) {    if (! numSamplesEqual(A,p.numQuad,elements->numElements) ) {
119          Finley_ErrorCode=TYPE_ERROR;          sprintf(error_msg,"__FILE__: sample points of coefficient A don't match (%d,%d)",p.numQuad,elements->numElements);
120          sprintf(Finley_ErrorMsg,"sample points of coefficient A don't match (%d,%d)",p.numQuad,elements->numElements);          Finley_setError(TYPE_ERROR,error_msg);
121    }    }
122    
123    if (! numSamplesEqual(B,p.numQuad,elements->numElements) ) {    if (! numSamplesEqual(B,p.numQuad,elements->numElements) ) {
124          Finley_ErrorCode=TYPE_ERROR;          sprintf(error_msg,"__FILE__: sample points of coefficient B don't match (%d,%d)",p.numQuad,elements->numElements);
125          sprintf(Finley_ErrorMsg,"sample points of coefficient B don't match (%d,%d)",p.numQuad,elements->numElements);          Finley_setError(TYPE_ERROR,error_msg);
126    }    }
127    
128    if (! numSamplesEqual(C,p.numQuad,elements->numElements) ) {    if (! numSamplesEqual(C,p.numQuad,elements->numElements) ) {
129          Finley_ErrorCode=TYPE_ERROR;          sprintf(error_msg,"__FILE__: sample points of coefficient C don't match (%d,%d)",p.numQuad,elements->numElements);
130          sprintf(Finley_ErrorMsg,"sample points of coefficient C don't match (%d,%d)",p.numQuad,elements->numElements);          Finley_setError(TYPE_ERROR,error_msg);
131    }    }
132    
133    if (! numSamplesEqual(D,p.numQuad,elements->numElements) ) {    if (! numSamplesEqual(D,p.numQuad,elements->numElements) ) {
134          Finley_ErrorCode=TYPE_ERROR;          sprintf(error_msg,"__FILE__: sample points of coefficient D don't match (%d,%d)",p.numQuad,elements->numElements);
135          sprintf(Finley_ErrorMsg,"sample points of coefficient D don't match (%d,%d)",p.numQuad,elements->numElements);          Finley_setError(TYPE_ERROR,error_msg);
136    }    }
137    
138    if (! numSamplesEqual(X,p.numQuad,elements->numElements) ) {    if (! numSamplesEqual(X,p.numQuad,elements->numElements) ) {
139          Finley_ErrorCode=TYPE_ERROR;          sprintf(error_msg,"__FILE__: sample points of coefficient X don't match (%d,%d)",p.numQuad,elements->numElements);
140          sprintf(Finley_ErrorMsg,"sample points of coefficient X don't match (%d,%d)",p.numQuad,elements->numElements);          Finley_setError(TYPE_ERROR,error_msg);
141    }    }
142    
143    if (! numSamplesEqual(Y,p.numQuad,elements->numElements) ) {    if (! numSamplesEqual(Y,p.numQuad,elements->numElements) ) {
144          Finley_ErrorCode=TYPE_ERROR;          sprintf(error_msg,"__FILE__: sample points of coefficient Y don't match (%d,%d)",p.numQuad,elements->numElements);
145          sprintf(Finley_ErrorMsg,"sample points of coefficient Y don't match (%d,%d)",p.numQuad,elements->numElements);          Finley_setError(TYPE_ERROR,error_msg);
146    }    }
147    
148    /*  check the dimensions: */    /*  check the dimensions: */
# Line 150  void Finley_Assemble_PDE(Finley_NodeFile Line 152  void Finley_Assemble_PDE(Finley_NodeFile
152        dimensions[0]=p.numDim;        dimensions[0]=p.numDim;
153        dimensions[1]=p.numDim;        dimensions[1]=p.numDim;
154        if (!isDataPointShapeEqual(A,2,dimensions)) {        if (!isDataPointShapeEqual(A,2,dimensions)) {
155            Finley_ErrorCode=TYPE_ERROR;            sprintf(error_msg,"__FILE__: coefficient A: illegal shape, expected shape (%d,%d)",dimensions[0],dimensions[1]);
156            sprintf(Finley_ErrorMsg,"coefficient A: illegal shape, expected shape (%d,%d)",dimensions[0],dimensions[1]);            Finley_setError(TYPE_ERROR,error_msg);
157        }        }
158      }      }
159      if (!isEmpty(B)) {      if (!isEmpty(B)) {
160         dimensions[0]=p.numDim;         dimensions[0]=p.numDim;
161         if (!isDataPointShapeEqual(B,1,dimensions)) {         if (!isDataPointShapeEqual(B,1,dimensions)) {
162            Finley_ErrorCode=TYPE_ERROR;            sprintf(error_msg,"__FILE__: coefficient B: illegal shape (%d,)",dimensions[0]);
163            sprintf(Finley_ErrorMsg,"coefficient B: illegal shape (%d,)",dimensions[0]);            Finley_setError(TYPE_ERROR,error_msg);
164         }         }
165      }      }
166      if (!isEmpty(C)) {      if (!isEmpty(C)) {
167         dimensions[0]=p.numDim;         dimensions[0]=p.numDim;
168         if (!isDataPointShapeEqual(C,1,dimensions)) {         if (!isDataPointShapeEqual(C,1,dimensions)) {
169            Finley_ErrorCode=TYPE_ERROR;            sprintf(error_msg,"__FILE__: coefficient C, expected shape (%d,)",dimensions[0]);
170            sprintf(Finley_ErrorMsg,"coefficient C, expected shape (%d,)",dimensions[0]);            Finley_setError(TYPE_ERROR,error_msg);
171         }         }
172      }      }
173      if (!isEmpty(D)) {      if (!isEmpty(D)) {
174         if (!isDataPointShapeEqual(D,0,dimensions)) {         if (!isDataPointShapeEqual(D,0,dimensions)) {
175            Finley_ErrorCode=TYPE_ERROR;            Finley_setError(TYPE_ERROR,"__FILE__: coefficient D, rank 0 expected.");
           sprintf(Finley_ErrorMsg,"coefficient D, rank 0 expected.");  
176         }         }
177      }      }
178      if (!isEmpty(X)) {      if (!isEmpty(X)) {
179         dimensions[0]=p.numDim;         dimensions[0]=p.numDim;
180         if (!isDataPointShapeEqual(X,1,dimensions)) {         if (!isDataPointShapeEqual(X,1,dimensions)) {
181            Finley_ErrorCode=TYPE_ERROR;            sprintf(error_msg,"__FILE__: coefficient X, expected shape (%d,",dimensions[0]);
182            sprintf(Finley_ErrorMsg,"coefficient X, expected shape (%d,",dimensions[0]);            Finley_setError(TYPE_ERROR,error_msg);
183         }         }
184      }      }
185      if (!isEmpty(Y)) {      if (!isEmpty(Y)) {
186         if (!isDataPointShapeEqual(Y,0,dimensions)) {         if (!isDataPointShapeEqual(Y,0,dimensions)) {
187            Finley_ErrorCode=TYPE_ERROR;            Finley_setError(TYPE_ERROR,"__FILE__: coefficient Y, rank 0 expected.");
           sprintf(Finley_ErrorMsg,"coefficient Y, rank 0 expected.");  
188         }         }
189      }      }
190    } else {    } else {
# Line 194  void Finley_Assemble_PDE(Finley_NodeFile Line 194  void Finley_Assemble_PDE(Finley_NodeFile
194        dimensions[2]=p.numComp;        dimensions[2]=p.numComp;
195        dimensions[3]=p.numDim;        dimensions[3]=p.numDim;
196        if (!isDataPointShapeEqual(A,4,dimensions)) {        if (!isDataPointShapeEqual(A,4,dimensions)) {
197            Finley_ErrorCode=TYPE_ERROR;            sprintf(error_msg,"__FILE__: coefficient A, expected shape (%d,%d,%d,%d)",dimensions[0],dimensions[1],dimensions[2],dimensions[3]);
198            sprintf(Finley_ErrorMsg,"coefficient A, expected shape (%d,%d,%d,%d)",dimensions[0],dimensions[1],dimensions[2],dimensions[3]);            Finley_setError(TYPE_ERROR,error_msg);
199        }        }
200      }      }
201      if (!isEmpty(B)) {      if (!isEmpty(B)) {
# Line 203  void Finley_Assemble_PDE(Finley_NodeFile Line 203  void Finley_Assemble_PDE(Finley_NodeFile
203        dimensions[1]=p.numDim;        dimensions[1]=p.numDim;
204        dimensions[2]=p.numComp;        dimensions[2]=p.numComp;
205        if (!isDataPointShapeEqual(B,3,dimensions)) {        if (!isDataPointShapeEqual(B,3,dimensions)) {
206            Finley_ErrorCode=TYPE_ERROR;            sprintf(error_msg,"__FILE__: coefficient B, expected shape (%d,%d,%d)",dimensions[0],dimensions[1],dimensions[2]);
207            sprintf(Finley_ErrorMsg,"coefficient B, expected shape (%d,%d,%d)",dimensions[0],dimensions[1],dimensions[2]);            Finley_setError(TYPE_ERROR,error_msg);
208        }        }
209      }      }
210      if (!isEmpty(C)) {      if (!isEmpty(C)) {
# Line 212  void Finley_Assemble_PDE(Finley_NodeFile Line 212  void Finley_Assemble_PDE(Finley_NodeFile
212        dimensions[1]=p.numComp;        dimensions[1]=p.numComp;
213        dimensions[2]=p.numDim;        dimensions[2]=p.numDim;
214        if (!isDataPointShapeEqual(C,3,dimensions)) {        if (!isDataPointShapeEqual(C,3,dimensions)) {
215             Finley_ErrorCode=TYPE_ERROR;            sprintf(error_msg,"__FILE__: coefficient C, expected shape (%d,%d,%d)",dimensions[0],dimensions[1],dimensions[2]);
216            sprintf(Finley_ErrorMsg,"coefficient C, expected shape (%d,%d,%d)",dimensions[0],dimensions[1],dimensions[2]);            Finley_setError(TYPE_ERROR,error_msg);
217        }        }
218      }      }
219      if (!isEmpty(D)) {      if (!isEmpty(D)) {
220        dimensions[0]=p.numEqu;        dimensions[0]=p.numEqu;
221        dimensions[1]=p.numComp;        dimensions[1]=p.numComp;
222        if (!isDataPointShapeEqual(D,2,dimensions)) {        if (!isDataPointShapeEqual(D,2,dimensions)) {
223            Finley_ErrorCode=TYPE_ERROR;            sprintf(error_msg,"__FILE__: coefficient D, expected shape (%d,%d)",dimensions[0],dimensions[1]);
224            sprintf(Finley_ErrorMsg,"coefficient D, expected shape (%d,%d)",dimensions[0],dimensions[1]);            Finley_setError(TYPE_ERROR,error_msg);
225        }        }
226      }      }
227      if (!isEmpty(X)) {      if (!isEmpty(X)) {
228        dimensions[0]=p.numEqu;        dimensions[0]=p.numEqu;
229        dimensions[1]=p.numDim;        dimensions[1]=p.numDim;
230        if (!isDataPointShapeEqual(X,2,dimensions)) {        if (!isDataPointShapeEqual(X,2,dimensions)) {
231             Finley_ErrorCode=TYPE_ERROR;            sprintf(error_msg,"__FILE__: coefficient X, expected shape (%d,%d)",dimensions[0],dimensions[1]);
232            sprintf(Finley_ErrorMsg,"coefficient X, expected shape (%d,%d)",dimensions[0],dimensions[1]);            Finley_setError(TYPE_ERROR,error_msg);
233        }        }
234      }      }
235      if (!isEmpty(Y)) {      if (!isEmpty(Y)) {
236        dimensions[0]=p.numEqu;        dimensions[0]=p.numEqu;
237        if (!isDataPointShapeEqual(Y,1,dimensions)) {        if (!isDataPointShapeEqual(Y,1,dimensions)) {
238             Finley_ErrorCode=TYPE_ERROR;            sprintf(error_msg,"__FILE__: coefficient Y, expected shape (%d,)",dimensions[0]);
239            sprintf(Finley_ErrorMsg,"coefficient Y, expected shape (%d,)",dimensions[0]);            Finley_setError(TYPE_ERROR,error_msg);
240        }        }
241      }      }
242    }    }
243    
244    if (Finley_ErrorCode==NO_ERROR) {    if (Finley_noError()) {
245       time0=Finley_timer();       time0=Finley_timer();
246       #pragma omp parallel private(index_col,index_row,EM_S,EM_F,V,dVdv,dSdV,Vol,dvdV,color,q) \       #pragma omp parallel private(index_col,index_row,EM_S,EM_F,V,dVdv,dSdV,Vol,dvdV,color,q) \
247              firstprivate(elements,nodes,S,F,A,B,C,D,X,Y)              firstprivate(elements,nodes,S,F,A,B,C,D,X,Y)
# Line 304  void Finley_Assemble_PDE(Finley_NodeFile Line 304  void Finley_Assemble_PDE(Finley_NodeFile
304                                                                       getSampleData(D,e),isExpanded(D));                                                                       getSampleData(D,e),isExpanded(D));
305                         }                         }
306                         for (q=0;q<p.NN_col;q++) index_col[q]=p.label_col[elements->Nodes[INDEX2(p.col_node[q],e,p.NN)]];                         for (q=0;q<p.NN_col;q++) index_col[q]=p.label_col[elements->Nodes[INDEX2(p.col_node[q],e,p.NN)]];
307                         Finley_SystemMatrix_add(S,p.NN_row,index_row,p.numEqu,p.NN_col,index_col,p.numComp,EM_S);                         Finley_Assemble_addToSystemMatrix(S,p.NN_row,index_row,p.numEqu,p.NN_col,index_col,p.numComp,EM_S);
308                       }                       }
309                       if (!isEmpty(F)) {                       if (!isEmpty(F)) {
310                         for (q=0;q<p.NN_row*p.numEqu;q++) EM_F[q]=0;                         for (q=0;q<p.NN_row*p.numEqu;q++) EM_F[q]=0;
# Line 344  void Finley_Assemble_PDE(Finley_NodeFile Line 344  void Finley_Assemble_PDE(Finley_NodeFile
344  }  }
345  /*  /*
346   * $Log$   * $Log$
347     * Revision 1.8  2005/09/15 03:44:21  jgs
348     * Merge of development branch dev-02 back to main trunk on 2005-09-15
349     *
350   * Revision 1.7  2005/09/01 03:31:35  jgs   * Revision 1.7  2005/09/01 03:31:35  jgs
351   * Merge of development branch dev-02 back to main trunk on 2005-09-01   * Merge of development branch dev-02 back to main trunk on 2005-09-01
352   *   *
353   * Revision 1.6  2005/08/12 01:45:42  jgs   * Revision 1.6  2005/08/12 01:45:42  jgs
354   * erge of development branch dev-02 back to main trunk on 2005-08-12   * erge of development branch dev-02 back to main trunk on 2005-08-12
355   *   *
356     * Revision 1.5.2.3  2005/09/07 06:26:17  gross
357     * the solver from finley are put into the standalone package paso now
358     *
359   * Revision 1.5.2.2  2005/08/24 02:02:18  gross   * Revision 1.5.2.2  2005/08/24 02:02:18  gross
360   * timing output switched off. solver output can be swiched through getSolution(verbose=True) now.   * timing output switched off. solver output can be swiched through getSolution(verbose=True) now.
361   *   *

Legend:
Removed from v.149  
changed lines
  Added in v.155

  ViewVC Help
Powered by ViewVC 1.1.26