Changeset 5854
- Timestamp:
- 09/16/10 17:38:47 (15 years ago)
- Location:
- issm/trunk/src/c/objects/Elements
- Files:
-
- 4 edited
Legend:
- Unmodified
- Added
- Removed
-
issm/trunk/src/c/objects/Elements/Penta.cpp
r5853 r5854 687 687 void Penta::CreateKMatrix(Mat Kgg, Mat Kff, Mat Kfs){ 688 688 689 /*retri ve parameters: */689 /*retrieve parameters: */ 690 690 int analysis_type; 691 ElementMatrix* Ke=NULL; 691 692 parameters->FindParam(&analysis_type,AnalysisTypeEnum); 692 693 … … 722 723 break; 723 724 case MeltingAnalysisEnum: 724 CreateKMatrixMelting( Kgg); 725 Ke=CreateKMatrixMelting(); 726 if(Ke) Ke->AddToGlobal(Kgg,NULL,NULL); 727 delete Ke; 725 728 break; 726 729 default: … … 2321 2324 bedrock, in which case we spawn a tria element using the 3 first grids, and use it to build 2322 2325 the stiffness matrix. */ 2323 if (!IsOnBed() ) return NULL;2326 if (!IsOnBed() || IsOnWater()) return NULL; 2324 2327 2325 2328 /*Depth Averaging B*/ … … 2821 2824 /*}}}*/ 2822 2825 /*FUNCTION Penta::CreateKMatrixMelting {{{1*/ 2823 void Penta::CreateKMatrixMelting(Mat Kgg){ 2824 2825 Tria* tria=NULL; 2826 2827 /*If on water, skip: */ 2828 if(IsOnWater())return; 2829 2830 if (!IsOnBed()){ 2831 return; 2832 } 2833 else{ 2834 2835 tria=(Tria*)SpawnTria(0,1,2); //grids 0, 1 and 2 make the new tria. 2836 tria->CreateKMatrixMelting(Kgg); 2837 delete tria->matice; delete tria; 2838 return; 2839 } 2826 ElementMatrix* Penta::CreateKMatrixMelting(void){ 2827 2828 if (!IsOnBed() || IsOnWater()) return NULL; 2829 2830 Tria* tria=(Tria*)SpawnTria(0,1,2); //grids 0, 1 and 2 make the new tria. 2831 ElementMatrix* Ke=tria->CreateKMatrixMelting(); 2832 2833 delete tria->matice; delete tria; 2834 return Ke; 2840 2835 } 2841 2836 /*}}}*/ -
issm/trunk/src/c/objects/Elements/Penta.h
r5853 r5854 140 140 ElementMatrix* CreateKMatrixDiagnosticStokesFriction(void); 141 141 void CreateKMatrixDiagnosticVert( Mat Kgg); 142 void CreateKMatrixMelting(Mat Kggg);142 ElementMatrix* CreateKMatrixMelting(void); 143 143 void CreateKMatrixPrognostic(Mat Kggg); 144 144 void CreateKMatrixSlope(Mat Kggg); -
issm/trunk/src/c/objects/Elements/Tria.cpp
r5852 r5854 3154 3154 /*}}}*/ 3155 3155 /*FUNCTION Tria::CreateKMatrixMelting {{{1*/ 3156 void Tria::CreateKMatrixMelting(Mat Kgg){ 3157 3158 /*indexing: */ 3159 int i,j; 3156 ElementMatrix* Tria::CreateKMatrixMelting(void){ 3160 3157 3161 3158 const int numdof=NUMVERTICES*NDOF1; 3162 int* doflist=NULL; 3163 3164 /*Grid data: */ 3159 int i,j,ig; 3165 3160 double xyz_list[NUMVERTICES][3]; 3166 3167 /*Material constants */3168 3161 double heatcapacity,latentheat; 3169 3170 /* gaussian points: */3171 int ig;3172 3162 GaussTria *gauss=NULL; 3173 3174 /*matrices: */3175 3163 double Jdet; 3176 3164 double D_scalar; … … 3180 3168 double Ke_gaussian[numdof][numdof]={0.0}; 3181 3169 3182 /*Recover constants of ice */ 3170 /*Initialize Element matrix and return if necessary*/ 3171 if(IsOnWater()) return NULL; 3172 ElementMatrix* Ke=this->NewElementMatrix(NoneApproximationEnum); 3173 3174 /*Retrieve all inputs and parameters*/ 3175 GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES); 3183 3176 latentheat=matpar->GetLatentHeat(); 3184 3177 heatcapacity=matpar->GetHeatCapacity(); 3185 3186 /* Get node coordinates and dof list: */3187 GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES);3188 GetDofList(&doflist,NoneApproximationEnum,GsetEnum);3189 3190 3178 3191 3179 /* Start looping on the number of gauss (nodes on the bedrock) */ … … 3195 3183 gauss->GaussPoint(ig); 3196 3184 3185 GetL(&L[0], &xyz_list[0][0], gauss,NDOF1); 3197 3186 GetJacobianDeterminant2d(&Jdet, &xyz_list[0][0], gauss); 3198 3199 /*Get L matrix : */3200 GetL(&L[0], &xyz_list[0][0], gauss,NDOF1);3201 3202 /*Calculate DL on gauss point */3203 3187 D_scalar=latentheat/heatcapacity*gauss->weight*Jdet; 3204 3188 3205 /* Do the triple product tL*D*L: */3206 3189 MatrixMultiply(&L[0],numdof,1,0,&D_scalar,1,1,0,&tLD[0],0); 3207 3190 MatrixMultiply(&tLD[0],numdof,1,0,&L[0],1,numdof,0,&Ke_gaussian[0][0],0); 3208 3191 3209 for(i=0;i<NUMVERTICES;i++){ 3210 for(j=0;j<NUMVERTICES;j++){ 3211 K_terms[i][j]+=Ke_gaussian[i][j]; 3212 } 3213 } 3214 } 3215 3216 /*Add Ke_gg to global matrix Kgg: */ 3217 MatSetValues(Kgg,numdof,doflist,numdof,doflist,(const double*)K_terms,ADD_VALUES); 3192 for(i=0;i<numdof;i++) for(j=0;j<numdof;j++) Ke->values[i*numdof+j]+=Ke_gaussian[i][j]; 3193 } 3218 3194 3219 3195 /*Clean up and return*/ 3220 3196 delete gauss; 3221 xfree((void**)&doflist);3197 return Ke; 3222 3198 } 3223 3199 /*}}}*/ -
issm/trunk/src/c/objects/Elements/Tria.h
r5852 r5854 133 133 ElementMatrix* CreateKMatrixDiagnosticHutter(void); 134 134 void CreateKMatrixDiagnosticSurfaceVert(Mat Kgg); 135 void CreateKMatrixMelting(Mat Kgg);135 ElementMatrix* CreateKMatrixMelting(void); 136 136 ElementMatrix* CreateKMatrixPrognostic(void); 137 137 ElementMatrix* CreateKMatrixPrognostic_CG(void);
Note:
See TracChangeset
for help on using the changeset viewer.