Changeset 5878


Ignore:
Timestamp:
09/17/10 17:29:07 (15 years ago)
Author:
seroussi
Message:

cleanup in tria

File:
1 edited

Legend:

Unmodified
Added
Removed

  • issm/trunk/src/c/objects/Elements/Tria.cpp

    r5877 r5878  
    30983098ElementMatrix* Tria::CreateKMatrixMelting(void){
    30993099
     3100        /*Constants*/
    31003101        const int  numdof=NUMVERTICES*NDOF1;
    3101         int i,j,ig;
     3102
     3103        /*Intermediaries */
     3104        int        i,j,ig;
     3105        double     heatcapacity,latentheat;
     3106        double     Jdet,D_scalar;
    31023107        double     xyz_list[NUMVERTICES][3];
    3103         double     heatcapacity,latentheat;
    3104         GaussTria *gauss=NULL;
    3105         double     Jdet;
    3106         double     D_scalar;
    3107         double     K_terms[numdof][numdof]={0.0};
    31083108        double     L[3];
    31093109        double     tLD[3];
    31103110        double     Ke_gaussian[numdof][numdof]={0.0};
     3111        GaussTria *gauss=NULL;
    31113112
    31123113        /*Initialize Element matrix and return if necessary*/
    … …  
    31273128                GetL(&L[0], &xyz_list[0][0], gauss,NDOF1);
    31283129                GetJacobianDeterminant2d(&Jdet, &xyz_list[0][0], gauss);
     3130
    31293131                D_scalar=latentheat/heatcapacity*gauss->weight*Jdet;
    31303132
    … …  
    31573159ElementMatrix* Tria::CreateKMatrixPrognostic_CG(void){
    31583160
    3159         /* local declarations */
     3161        /*Constants*/
    31603162        const int    numdof=NDOF1*NUMVERTICES;
    3161         int             i,j;
    3162         int     ig;
    3163         double       xyz_list[NUMVERTICES][3];
    3164         int          numberofdofspernode=1;
     3163
     3164        /*Intermediaries */
     3165        bool       artdiff;
     3166        int        i,j,ig,dim;
     3167        double     Jdettria,DL_scalar,dt;
     3168        double     vx,vy,dvxdx,dvydy;
     3169        double     dvx[2],dvy[2];
     3170        double     v_gauss[2]={0.0};
     3171        double     xyz_list[NUMVERTICES][3];
     3172        double     L[NUMVERTICES];
     3173        double     B[2][NUMVERTICES];
     3174        double     Bprime[2][NUMVERTICES];
     3175        double     K[2][2]                        ={0.0};
     3176        double     KDL[2][2]                      ={0.0};
     3177        double     DL[2][2]                        ={0.0};
     3178        double     DLprime[2][2]                   ={0.0};
     3179        double     Ke_gg_gaussian[numdof][numdof]  ={0.0};
     3180        double     Ke_gg_thickness1[numdof][numdof]={0.0};
     3181        double     Ke_gg_thickness2[numdof][numdof]={0.0};
    31653182        GaussTria *gauss=NULL;
    3166         double L[NUMVERTICES];
    3167         double B[2][NUMVERTICES];
    3168         double Bprime[2][NUMVERTICES];
    3169         double DL[2][2]={0.0};
    3170         double DLprime[2][2]={0.0};
    3171         double DL_scalar;
    3172         double Ke_gg[numdof][numdof]={0.0};
    3173         double Ke_gg_gaussian[numdof][numdof]={0.0};
    3174         double Ke_gg_thickness1[numdof][numdof]={0.0};
    3175         double Ke_gg_thickness2[numdof][numdof]={0.0};
    3176         double Jdettria;
    3177         double  dvx[2];
    3178         double  dvy[2];
    3179         double  vx,vy;
    3180         double  dvxdx,dvydy;
    3181         double  v_gauss[2]={0.0};
    3182         double  K[2][2]={0.0};
    3183         double  KDL[2][2]={0.0};
    3184         int     dim;
    3185         double dt;
    3186         bool artdiff;
    31873183
    31883184        /*Initialize Element matrix and return if necessary*/
    … …  
    32083204        //Create Artificial diffusivity once for all if requested
    32093205        if(artdiff){
    3210                 //Get the Jacobian determinant
    32113206                gauss=new GaussTria();
    32123207                gauss->GaussCenter();
    … …  
    32143209                delete gauss;
    32153210
    3216                 //Build K matrix (artificial diffusivity matrix)
    32173211                vxaverage_input->GetParameterAverage(&v_gauss[0]);
    32183212                vyaverage_input->GetParameterAverage(&v_gauss[1]);
    … …  
    32283222                gauss->GaussPoint(ig);
    32293223
    3230                 GetL(&L[0], &xyz_list[0][0], gauss,numberofdofspernode);
    32313224                GetJacobianDeterminant2d(&Jdettria, &xyz_list[0][0],gauss);
     3225                GetL(&L[0], &xyz_list[0][0], gauss,NDOF1);
     3226
     3227                vxaverage_input->GetParameterValue(&vx,gauss);
     3228                vyaverage_input->GetParameterValue(&vy,gauss);
     3229                vxaverage_input->GetParameterDerivativeValue(&dvx[0],&xyz_list[0][0],gauss);
     3230                vyaverage_input->GetParameterDerivativeValue(&dvy[0],&xyz_list[0][0],gauss);
    32323231
    32333232                DL_scalar=gauss->weight*Jdettria;
    … …  
    32413240                GetBprimePrognostic(&Bprime[0][0], &xyz_list[0][0], gauss);
    32423241
    3243                 vxaverage_input->GetParameterValue(&vx,gauss);
    3244                 vyaverage_input->GetParameterValue(&vy,gauss);
    3245                 vxaverage_input->GetParameterDerivativeValue(&dvx[0],&xyz_list[0][0],gauss);
    3246                 vyaverage_input->GetParameterDerivativeValue(&dvy[0],&xyz_list[0][0],gauss);
    32473242                dvxdx=dvx[0];
    32483243                dvydy=dvy[1];
    … …  
    32893284ElementMatrix* Tria::CreateKMatrixPrognostic_DG(void){
    32903285
    3291         /* local declarations */
    3292         int             i,j;
     3286        /*Constants*/
    32933287        const int    numdof=NDOF1*NUMVERTICES;
    3294         double       xyz_list[NUMVERTICES][3];
    3295         int*         doflist=NULL;
    3296         int          numberofdofspernode=1;
    3297         int     ig;
    3298         GaussTria *gauss=NULL;
    3299         double L[NUMVERTICES];
    3300         double B[2][NUMVERTICES];
    3301         double Bprime[2][NUMVERTICES];
    3302         double DL[2][2]={0.0};
    3303         double DLprime[2][2]={0.0};
    3304         double DL_scalar;
    3305         double Ke_gg[numdof][numdof]={0.0};
    3306         double Ke_gg1[numdof][numdof]={0.0};
    3307         double Ke_gg2[numdof][numdof]={0.0};
    3308         double Jdettria;
    3309         double  vx,vy;
    3310         int     dim;
    3311         double dt;
     3288
     3289        /*Intermediaries */
     3290        int        i,j,ig,dim;
     3291        int*       doflist=NULL;
     3292        double     xyz_list[NUMVERTICES][3];
     3293        double     Jdettria,dt,vx,vy;
     3294        double     L[NUMVERTICES];
     3295        double     B[2][NUMVERTICES];
     3296        double     Bprime[2][NUMVERTICES];
     3297        double     DL[2][2]={0.0};
     3298        double     DLprime[2][2]={0.0};
     3299        double     DL_scalar;
     3300        double     Ke_gg1[numdof][numdof]={0.0};
     3301        double     Ke_gg2[numdof][numdof]={0.0};
     3302        GaussTria  *gauss=NULL;
    33123303
    33133304        /*Initialize Element matrix and return if necessary*/
    … …  
    33363327                gauss->GaussPoint(ig);
    33373328
    3338 
    3339                 GetL(&L[0], &xyz_list[0][0], gauss,numberofdofspernode);
     3329                vxaverage_input->GetParameterValue(&vx,gauss);
     3330                vyaverage_input->GetParameterValue(&vy,gauss);
     3331
    33403332                GetJacobianDeterminant2d(&Jdettria, &xyz_list[0][0],gauss);
     3333                GetL(&L[0], &xyz_list[0][0], gauss,NDOF1);
    33413334
    33423335                DL_scalar=gauss->weight*Jdettria;
    … …  
    33473340                                        &Ke_gg1[0][0],0);
    33483341
    3349                 /*Get B  and B prime matrix: */
    33503342                /*WARNING: B and Bprime are inverted compared to usual prognostic!!!!*/
    33513343                GetBPrognostic(&Bprime[0][0], &xyz_list[0][0], gauss);
    33523344                GetBprimePrognostic(&B[0][0], &xyz_list[0][0], gauss);
    33533345
    3354                 vxaverage_input->GetParameterValue(&vx,gauss);
    3355                 vyaverage_input->GetParameterValue(&vy,gauss);
    33563346                DL_scalar=-dt*gauss->weight*Jdettria;
    33573347
    … …  
    36323622void  Tria::CreatePVectorBalancedvelocities(Vec pg){
    36333623
    3634         /* local declarations */
    3635         int             i,j;
    3636 
    3637         /* node data: */
     3624        /*Constants*/
    36383625        const int    numdof=NDOF1*NUMVERTICES;
    3639         double       xyz_list[NUMVERTICES][3];
    3640         int*         doflist=NULL;
    3641         int          numberofdofspernode=1;
    3642 
    3643         /* gaussian points: */
    3644         int     ig;
     3626
     3627        /*Intermediaries */
     3628        int        i,j,ig;
     3629        int*       doflist=NULL;
     3630        double     xyz_list[NUMVERTICES][3];
     3631        double     Jdettria,accumulation_g,melting_g;
     3632        double     L[NUMVERTICES];
     3633        double     pe_g[NUMVERTICES]={0.0};
    36453634        GaussTria* gauss=NULL;
    36463635
    3647         /* matrix */
    3648         double pe_g[NUMVERTICES]={0.0};
    3649         double L[NUMVERTICES];
    3650         double Jdettria;
    3651 
    3652         /*input parameters for structural analysis (diagnostic): */
    3653         double  accumulation_g;
    3654         double  melting_g;
    3655 
    3656         /* Get node coordinates and dof list: */
    36573636        GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES);
    36583637        GetDofList(&doflist,NoneApproximationEnum,GsetEnum);
    … …  
    36683647                gauss->GaussPoint(ig);
    36693648
    3670                 /* Get Jacobian determinant: */
    36713649                GetJacobianDeterminant2d(&Jdettria, &xyz_list[0][0],gauss);
    3672 
    3673                 /*Get L matrix: */
    3674                 GetL(&L[0], &xyz_list[0][0], gauss,numberofdofspernode);
    3675 
    3676                 /* Get accumulation, melting at gauss point */
     3650                GetL(&L[0], &xyz_list[0][0], gauss,NDOF1);
     3651
    36773652                accumulation_input->GetParameterValue(&accumulation_g,gauss);
    36783653                melting_input->GetParameterValue(&melting_g,gauss);
    36793654
    3680                 /* Add value into pe_g: */
    36813655                for( i=0; i<numdof; i++) pe_g[i]+=Jdettria*gauss->weight*(accumulation_g-melting_g)*L[i];
    36823656
    36833657        }
    36843658
    3685         /*Add pe_g to global matrix Kgg: */
    36863659        VecSetValues(pg,numdof,doflist,(const double*)pe_g,ADD_VALUES);
    36873660
Note: See TracChangeset for help on using the changeset viewer.