Changeset 10523


Ignore:
Timestamp:
11/08/11 10:07:23 (13 years ago)
Author:
Mathieu Morlighem
Message:

Extended change of Coordinate Systems for coupling (not done yet)

Location:
issm/trunk/src/c
Files:
10 edited

Legend:

Unmodified
Added
Removed

  • issm/trunk/src/c/EnumDefinitions/EnumDefinitions.h

    r10521 r10523  
    445445        NearestInterpEnum,
    446446        /*}}}*/
     447        /*Coordinate Systems{{{1*/
     448        XYEnum,
     449        XYZPEnum,
     450        /*}}}*/
    447451        /*Options{{{1*/
    448452        OptionEnum,

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

    r10517 r10523  
    57455745
    57465746        /*Transform Coordinate System*/
    5747         TransformStiffnessMatrixCoord(Ke,tria->nodes,NUMVERTICES2D,NDOF2);
     5747        TransformStiffnessMatrixCoord(Ke,tria->nodes,NUMVERTICES2D,XYEnum);
    57485748
    57495749        /*Clean up and return*/
    … …  
    58775877
    58785878        /*Transform Coordinate System*/
    5879         TransformStiffnessMatrixCoord(Ke,nodes,NUMVERTICES,NDOF2);
     5879        TransformStiffnessMatrixCoord(Ke,nodes,NUMVERTICES,XYEnum);
    58805880
    58815881        /*Clean up and return*/
    … …  
    59515951
    59525952        /*Transform Coordinate System*/
    5953         TransformStiffnessMatrixCoord(Ke,nodes,NUMVERTICES,NDOF2);
     5953        TransformStiffnessMatrixCoord(Ke,nodes,NUMVERTICES,XYEnum);
    59545954
    59555955        /*Clean up and return*/
    … …  
    60436043
    60446044        /*Transform Coordinate System*/
    6045         TransformStiffnessMatrixCoord(Ke,nodes,NUMVERTICES,NDOF4);
     6045        TransformStiffnessMatrixCoord(Ke,nodes,NUMVERTICES,XYZPEnum);
    60466046
    60476047        /*Clean up and return*/
    … …  
    61146114
    61156115        /*DO NOT Transform Coordinate System: this stiffness matrix is already expressed in tangential coordinates*/
    6116         //TransformStiffnessMatrixCoord(Ke,nodes,NUMVERTICES,NDOF4);
     6116        //TransformStiffnessMatrixCoord(Ke,nodes,NUMVERTICES,XYZPEnum);
    61176117       
    61186118        /*Clean up and return*/
    … …  
    67356735
    67366736        /*Transform coordinate system*/
    6737         TransformLoadVectorCoord(pe,nodes,NUMVERTICES,NDOF2);
     6737        TransformLoadVectorCoord(pe,nodes,NUMVERTICES,XYEnum);
    67386738
    67396739        /*Clean up and return*/
    … …  
    68296829
    68306830        /*Transform coordinate system*/
    6831         TransformLoadVectorCoord(pe,nodes,NUMVERTICES,NDOF4);
     6831        TransformLoadVectorCoord(pe,nodes,NUMVERTICES,XYZPEnum);
    68326832
    68336833        /*Clean up and return*/
    … …  
    68976897
    68986898        /*Transform coordinate system*/
    6899         TransformLoadVectorCoord(pe,nodes,NUMVERTICES,NDOF4);
     6899        TransformLoadVectorCoord(pe,nodes,NUMVERTICES,XYZPEnum);
    69006900
    69016901        /*Clean up and return*/
    … …  
    72567256
    72577257        /*Transform solution in Cartesian Space*/
    7258         TransformSolutionCoord(&values[0],nodes,NUMVERTICES2D,NDOF2); /*2D: only the first 3 nodes are taken*/
     7258        TransformSolutionCoord(&values[0],nodes,NUMVERTICES2D,XYEnum); /*2D: only the first 3 nodes are taken*/
    72597259
    72607260        /*Ok, we have vx and vy in values, fill in vx and vy arrays and extrude */
    … …  
    75157515
    75167516        /*Transform solution in Cartesian Space*/
    7517         TransformSolutionCoord(&values[0],nodes,NUMVERTICES,NDOF2);
     7517        TransformSolutionCoord(&values[0],nodes,NUMVERTICES,XYEnum);
    75187518
    75197519        /*Ok, we have vx and vy in values, fill in vx and vy arrays: */
    … …  
    78417841
    78427842        /*Transform solution in Cartesian Space*/
    7843         TransformSolutionCoord(&values[0],nodes,NUMVERTICES,NDOF4);
     7843        TransformSolutionCoord(&values[0],nodes,NUMVERTICES,XYZPEnum);
    78447844
    78457845        /*Ok, we have vx and vy in values, fill in all arrays: */

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

    r10449 r10523  
    27072707
    27082708        /*Transform Coordinate System*/
    2709         TransformStiffnessMatrixCoord(Ke,nodes,NUMVERTICES,NDOF2);
     2709        TransformStiffnessMatrixCoord(Ke,nodes,NUMVERTICES,XYEnum);
    27102710
    27112711        /*Clean up and return*/
    … …  
    27752775
    27762776        /*Transform Coordinate System*/
    2777         TransformStiffnessMatrixCoord(Ke,nodes,NUMVERTICES,NDOF2);
     2777        TransformStiffnessMatrixCoord(Ke,nodes,NUMVERTICES,XYEnum);
    27782778
    27792779        /*Clean up and return*/
    … …  
    28512851
    28522852        /*Transform coordinate system*/
    2853         TransformLoadVectorCoord(pe,nodes,NUMVERTICES,NDOF2);
     2853        TransformLoadVectorCoord(pe,nodes,NUMVERTICES,XYEnum);
    28542854
    28552855        /*Clean up and return*/
    … …  
    30093009
    30103010        /*Transform solution in Cartesian Space*/
    3011         TransformSolutionCoord(&values[0],nodes,NUMVERTICES,NDOF2);
     3011        TransformSolutionCoord(&values[0],nodes,NUMVERTICES,XYEnum);
    30123012
    30133013        /*Ok, we have vx and vy in values, fill in vx and vy arrays: */

  • issm/trunk/src/c/objects/Loads/Icefront.cpp

    r10407 r10523  
    530530
    531531        /*Transform load vector*/
    532         TransformLoadVectorCoord(pe,nodes,NUMVERTICESSEG,NDOF2);
     532        TransformLoadVectorCoord(pe,nodes,NUMVERTICESSEG,XYEnum);
    533533
    534534        /*Clean up and return*/
    … …  
    647647
    648648        /*Transform load vector*/
    649         TransformLoadVectorCoord(pe,nodes,NUMVERTICESQUA,NDOF2);
     649        TransformLoadVectorCoord(pe,nodes,NUMVERTICESQUA,XYEnum);
    650650
    651651        /*Clean up and return*/
    … …  
    725725
    726726        /*Transform load vector*/
    727         TransformLoadVectorCoord(pe,nodes,NUMVERTICESQUA,NDOF4);
     727        TransformLoadVectorCoord(pe,nodes,NUMVERTICESQUA,XYZPEnum);
    728728
    729729        /*Clean up and return*/

  • issm/trunk/src/c/objects/Loads/Pengrid.cpp

    r10407 r10523  
    560560
    561561        /*Transform Coordinate System*/
    562         TransformStiffnessMatrixCoord(Ke,&node,NUMVERTICES,NDOF4);
     562        TransformStiffnessMatrixCoord(Ke,&node,NUMVERTICES,XYZPEnum);
    563563
    564564        /*Clean up and return*/

  • issm/trunk/src/c/shared/Elements/CoordinateSystemTransform.cpp

    r10407 r10523  
    55#include <math.h>
    66
    7 void CoordinateSystemTransform(double** ptransform,Node** nodes,int numnodes,int dofspernode){
     7void CoordinateSystemTransform(double** ptransform,Node** nodes,int numnodes,int* cs_array){
    88
    9         int     i;
     9        int     i,counter;
     10        int     numdofs           = 0;
    1011        double  x_norm,y_norm;
    1112        double *transform         = NULL;
    … …  
    1314        double  coord_system[3][3];
    1415
    15         /*Some checks*/
     16        /*Some checks in debugging mode*/
    1617        _assert_(numnodes && nodes);
    17         if(dofspernode!=2 && dofspernode!=3 && dofspernode!=4) _error_("list of dofs per node supported: 2, 3 and 4");
     18
     19        /*Get total number of dofs*/
     20        for(i=0;i<numnodes;i++){
     21                switch(cs_array[i]){
     22                        case XYEnum:   numdofs+=2; break;
     23                        case XYZPEnum: numdofs+=4; break;
     24                        default: _error_("Coordinate system %s not supported yet",EnumToStringx(cs_array[i]));
     25                }
     26        }
    1827
    1928        /*Allocate and initialize transform matrix*/
    20         transform=(double*)xmalloc(numnodes*dofspernode*numnodes*dofspernode*sizeof(double));
    21         for(i=0;i<numnodes*dofspernode*numnodes*dofspernode;i++) transform[i]=0.0;
     29        transform=(double*)xmalloc(numdofs*numdofs*sizeof(double));
     30        for(i=0;i<numdofs*numdofs;i++) transform[i]=0.0;
    2231
    2332        /*Create transform matrix for all nodes (x,y for 2d and x,y,z for 3d). It is a block matrix
    … …  
    3039         * Where T1 is the transform matrix for node 1. It is a simple copy of the coordinate system
    3140         * associated to this node*/
     41        counter=0;
    3242        for(i=0;i<numnodes;i++){
    3343                nodes[i]->GetCoordinateSystem(&coord_system[0][0]);
    34                 if(dofspernode==2){
    35                         /*We remove the z component, we need to renormalize x and y*/
    36                         x_norm = sqrt( coord_system[0][0]*coord_system[0][0] + coord_system[1][0]*coord_system[1][0]);
    37                         y_norm = sqrt( coord_system[0][1]*coord_system[0][1] + coord_system[1][1]*coord_system[1][1]);
    38                         transform[(dofspernode*numnodes)*(i*dofspernode+0) + i*dofspernode+0] = coord_system[0][0]/x_norm;
    39                         transform[(dofspernode*numnodes)*(i*dofspernode+0) + i*dofspernode+1] = coord_system[0][1]/y_norm;
    40                         transform[(dofspernode*numnodes)*(i*dofspernode+1) + i*dofspernode+0] = coord_system[1][0]/x_norm;
    41                         transform[(dofspernode*numnodes)*(i*dofspernode+1) + i*dofspernode+1] = coord_system[1][1]/y_norm;
    42                 }
    43                 else if(dofspernode==3){
    44                         transform[(dofspernode*numnodes)*(i*dofspernode+0) + i*dofspernode+0] = coord_system[0][0];
    45                         transform[(dofspernode*numnodes)*(i*dofspernode+0) + i*dofspernode+1] = coord_system[0][1];
    46                         transform[(dofspernode*numnodes)*(i*dofspernode+0) + i*dofspernode+2] = coord_system[0][2];
    47                         transform[(dofspernode*numnodes)*(i*dofspernode+1) + i*dofspernode+0] = coord_system[1][0];
    48                         transform[(dofspernode*numnodes)*(i*dofspernode+1) + i*dofspernode+1] = coord_system[1][1];
    49                         transform[(dofspernode*numnodes)*(i*dofspernode+1) + i*dofspernode+2] = coord_system[1][2];
    50                         transform[(dofspernode*numnodes)*(i*dofspernode+2) + i*dofspernode+0] = coord_system[2][0];
    51                         transform[(dofspernode*numnodes)*(i*dofspernode+2) + i*dofspernode+1] = coord_system[2][1];
    52                         transform[(dofspernode*numnodes)*(i*dofspernode+2) + i*dofspernode+2] = coord_system[2][2];
    53                 }
    54                 else if(dofspernode==4){
    55                         /*Only the first 3 coordinates are changed (x,y,z), leave the others (P) unchanged*/
    56                         transform[(dofspernode*numnodes)*(i*dofspernode+0) + i*dofspernode+0] = coord_system[0][0];
    57                         transform[(dofspernode*numnodes)*(i*dofspernode+0) + i*dofspernode+1] = coord_system[0][1];
    58                         transform[(dofspernode*numnodes)*(i*dofspernode+0) + i*dofspernode+2] = coord_system[0][2];
    59                         transform[(dofspernode*numnodes)*(i*dofspernode+1) + i*dofspernode+0] = coord_system[1][0];
    60                         transform[(dofspernode*numnodes)*(i*dofspernode+1) + i*dofspernode+1] = coord_system[1][1];
    61                         transform[(dofspernode*numnodes)*(i*dofspernode+1) + i*dofspernode+2] = coord_system[1][2];
    62                         transform[(dofspernode*numnodes)*(i*dofspernode+2) + i*dofspernode+0] = coord_system[2][0];
    63                         transform[(dofspernode*numnodes)*(i*dofspernode+2) + i*dofspernode+1] = coord_system[2][1];
    64                         transform[(dofspernode*numnodes)*(i*dofspernode+2) + i*dofspernode+2] = coord_system[2][2];
    65                         transform[(dofspernode*numnodes)*(i*dofspernode+3) + i*dofspernode+3] = 1.0;
     44                switch(cs_array[i]){
     45                        case XYEnum:
     46                                /*We remove the z component, we need to renormalize x and y*/
     47                                x_norm = sqrt( coord_system[0][0]*coord_system[0][0] + coord_system[1][0]*coord_system[1][0]);
     48                                y_norm = sqrt( coord_system[0][1]*coord_system[0][1] + coord_system[1][1]*coord_system[1][1]);
     49                                transform[(numdofs)*(counter+0) + counter+0] = coord_system[0][0]/x_norm;
     50                                transform[(numdofs)*(counter+0) + counter+1] = coord_system[0][1]/y_norm;
     51                                transform[(numdofs)*(counter+1) + counter+0] = coord_system[1][0]/x_norm;
     52                                transform[(numdofs)*(counter+1) + counter+1] = coord_system[1][1]/y_norm;
     53                                counter+=2;
     54                                break;
     55                        case XYZPEnum:
     56                                /*Only the first 3 coordinates are changed (x,y,z), leave the others (P) unchanged*/
     57                                transform[(numdofs)*(counter+0) + counter+0] = coord_system[0][0];
     58                                transform[(numdofs)*(counter+0) + counter+1] = coord_system[0][1];
     59                                transform[(numdofs)*(counter+0) + counter+2] = coord_system[0][2];
     60                                transform[(numdofs)*(counter+1) + counter+0] = coord_system[1][0];
     61                                transform[(numdofs)*(counter+1) + counter+1] = coord_system[1][1];
     62                                transform[(numdofs)*(counter+1) + counter+2] = coord_system[1][2];
     63                                transform[(numdofs)*(counter+2) + counter+0] = coord_system[2][0];
     64                                transform[(numdofs)*(counter+2) + counter+1] = coord_system[2][1];
     65                                transform[(numdofs)*(counter+2) + counter+2] = coord_system[2][2];
     66                                transform[(numdofs)*(counter+3) + counter+3] = 1.0;
     67                                counter+=4;
     68                                break;
     69                        default:
     70                                _error_("Coordinate system %s not supported yet",EnumToStringx(cs_array[i]));
    6671                }
    6772        }

  • issm/trunk/src/c/shared/Elements/TransformLoadVectorCoord.cpp

    r10407 r10523  
    44#include "./elements.h"
    55
    6 void TransformLoadVectorCoord(ElementVector* pe,Node** nodes,int numnodes,int dofspernode){
     6void TransformLoadVectorCoord(ElementVector* pe,Node** nodes,int numnodes,int cs_enum){
     7        int* cs_array=NULL;
     8
     9        /*All nodes have the same Coordinate System*/
     10        cs_array=(int*)xmalloc(numnodes*sizeof(int));
     11        for(int i=0;i<numnodes;i++) cs_array[i]=cs_enum;
     12
     13        /*Call core*/
     14        TransformLoadVectorCoord(pe,nodes,numnodes,cs_array);
     15
     16        /*Clean-up*/
     17        xfree((void**)&cs_array);
     18}
     19
     20void TransformLoadVectorCoord(ElementVector* pe,Node** nodes,int numnodes,int* cs_array){
    721
    822        int     i,j;
     23        int     numdofs   = 0;
    924        double *transform = NULL;
    1025        double *values    = NULL;
     26
     27        /*Get total number of dofs*/
     28        for(i=0;i<numnodes;i++){
     29                switch(cs_array[i]){
     30                        case XYEnum:   numdofs+=2; break;
     31                        case XYZPEnum: numdofs+=4; break;
     32                        default: _error_("Coordinate system %s not supported yet",EnumToStringx(cs_array[i]));
     33                }
     34        }
    1135
    1236        /*Copy current load vector*/
    … …  
    1539
    1640        /*Get Coordinate Systems transform matrix*/
    17         CoordinateSystemTransform(&transform,nodes,numnodes,dofspernode);
     41        CoordinateSystemTransform(&transform,nodes,numnodes,cs_array);
    1842
    1943        /*Transform matrix: R^T*pe */
    20         MatrixMultiply(transform,numnodes*dofspernode,numnodes*dofspernode,1,
     44        MatrixMultiply(transform,numdofs,numdofs,1,
    2145                                values,pe->nrows,1,0,
    2246                                &pe->values[0],0);

  • issm/trunk/src/c/shared/Elements/TransformSolutionCoord.cpp

    r10407 r10523  
    44#include "./elements.h"
    55
    6 void TransformSolutionCoord(double* solution,Node** nodes,int numnodes,int dofspernode){
     6void TransformSolutionCoord(double* solution,Node** nodes,int numnodes,int cs_enum){
     7
     8        int* cs_array=NULL;
     9
     10        /*All nodes have the same Coordinate System*/
     11        cs_array=(int*)xmalloc(numnodes*sizeof(int));
     12        for(int i=0;i<numnodes;i++) cs_array[i]=cs_enum;
     13
     14        /*Call core*/
     15        TransformSolutionCoord(solution,nodes,numnodes,cs_array);
     16
     17        /*Clean-up*/
     18        xfree((void**)&cs_array);
     19}
     20
     21void TransformSolutionCoord(double* solution,Node** nodes,int numnodes,int* cs_array){
    722
    823        int     i,j;
     24        int     numdofs   = 0;
    925        double *transform = NULL;
    1026        double *values    = NULL;
    1127
     28        /*Get total number of dofs*/
     29        for(i=0;i<numnodes;i++){
     30                switch(cs_array[i]){
     31                        case XYEnum:   numdofs+=2; break;
     32                        case XYZPEnum: numdofs+=4; break;
     33                        default: _error_("Coordinate system %s not supported yet",EnumToStringx(cs_array[i]));
     34                }
     35        }
     36
    1237        /*Copy current solution vector*/
    13         values=(double*)xmalloc(numnodes*dofspernode*sizeof(double));
    14         for(i=0;i<numnodes*dofspernode;i++) values[i]=solution[i];
     38        values=(double*)xmalloc(numdofs*sizeof(double));
     39        for(i=0;i<numdofs;i++) values[i]=solution[i];
    1540
    1641        /*Get Coordinate Systems transform matrix*/
    17         CoordinateSystemTransform(&transform,nodes,numnodes,dofspernode);
     42        CoordinateSystemTransform(&transform,nodes,numnodes,cs_array);
    1843
    1944        /*Transform matrix: R*U */
    20         MatrixMultiply(transform,numnodes*dofspernode,numnodes*dofspernode,0,
    21                                 values,numnodes*dofspernode,1,0,
     45        MatrixMultiply(transform,numdofs,numdofs,0,
     46                                values,numdofs,1,0,
    2247                                &solution[0],0);
    2348

  • issm/trunk/src/c/shared/Elements/TransformStiffnessMatrixCoord.cpp

    r10407 r10523  
    44#include "./elements.h"
    55
    6 void   TransformStiffnessMatrixCoord(ElementMatrix* Ke,Node** nodes,int numnodes,int dofspernode){
     6void TransformStiffnessMatrixCoord(ElementMatrix* Ke,Node** nodes,int numnodes,int cs_enum){
     7
     8        int* cs_array=NULL;
     9
     10        /*All nodes have the same Coordinate System*/
     11        cs_array=(int*)xmalloc(numnodes*sizeof(int));
     12        for(int i=0;i<numnodes;i++) cs_array[i]=cs_enum;
     13
     14        /*Call core*/
     15        TransformStiffnessMatrixCoord(Ke,nodes,numnodes,cs_array);
     16
     17        /*Clean-up*/
     18        xfree((void**)&cs_array);
     19}
     20
     21void TransformStiffnessMatrixCoord(ElementMatrix* Ke,Node** nodes,int numnodes,int* cs_array){
    722
    823        int     i,j;
     24        int     numdofs   = 0;
    925        double *transform = NULL;
    1026        double *values    = NULL;
     27
     28        /*Get total number of dofs*/
     29        for(i=0;i<numnodes;i++){
     30                switch(cs_array[i]){
     31                        case XYEnum:   numdofs+=2; break;
     32                        case XYZPEnum: numdofs+=4; break;
     33                        default: _error_("Coordinate system %s not supported yet",EnumToStringx(cs_array[i]));
     34                }
     35        }
    1136
    1237        /*Copy current stiffness matrix*/
    … …  
    1540
    1641        /*Get Coordinate Systems transform matrix*/
    17         CoordinateSystemTransform(&transform,nodes,numnodes,dofspernode);
     42        CoordinateSystemTransform(&transform,nodes,numnodes,cs_array);
    1843
    1944        /*Transform matrix: R^T*Ke*R */
    20         TripleMultiply(transform,numnodes*dofspernode,numnodes*dofspernode,1,
     45        TripleMultiply(transform,numdofs,numdofs,1,
    2146                                values,Ke->nrows,Ke->ncols,0,
    22                                 transform,numnodes*dofspernode,numnodes*dofspernode,0,
     47                                transform,numdofs,numdofs,0,
    2348                                &Ke->values[0],0);
    2449

  • issm/trunk/src/c/shared/Elements/elements.h

    r10480 r10523  
    1717int*   GetLocalDofList( Node** nodes,int numnodes,int setenum,int approximation_enum);
    1818int*   GetGlobalDofList(Node** nodes,int numnodes,int setenum,int approximation_enum);
    19 void   CoordinateSystemTransform(double** ptransform,Node** nodes,int numnodes,int dofspernode);
     19void   CoordinateSystemTransform(double** ptransform,Node** nodes,int numnodes,int* cs_array);
    2020#ifdef _HAVE_DIAGNOSTIC_
    21 void   TransformStiffnessMatrixCoord(ElementMatrix* Ke,Node** nodes,int numnodes,int dofspernode);
    22 void   TransformLoadVectorCoord(ElementVector* pe,Node** nodes,int numnodes,int dofspernode);
    23 void   TransformSolutionCoord(double* solution,Node** nodes,int numnodes,int dofspernode);
     21void   TransformStiffnessMatrixCoord(ElementMatrix* Ke,Node** nodes,int numnodes,int cs_enum);
     22void   TransformStiffnessMatrixCoord(ElementMatrix* Ke,Node** nodes,int numnodes,int* cs_array);
     23void   TransformLoadVectorCoord(ElementVector* pe,Node** nodes,int numnodes,int cs_enum);
     24void   TransformLoadVectorCoord(ElementVector* pe,Node** nodes,int numnodes,int* cs_array);
     25void   TransformSolutionCoord(double* solution,Node** nodes,int numnodes,int cs_enum);
     26void   TransformSolutionCoord(double* solution,Node** nodes,int numnodes,int* cs_array);
    2427#endif
    2528
Note: See TracChangeset for help on using the changeset viewer.