Changeset 15635


Ignore:
Timestamp:
07/26/13 08:44:23 (12 years ago)
Author:
Mathieu Morlighem
Message:

NEW: added support for P2xP1 constraints

File:
1 edited

Legend:

Unmodified
Added
Removed

  • issm/trunk-jpl/src/c/modules/IoModelToConstraintsx/IoModelToConstraintsx.cpp

    r15630 r15635  
    3838        iomodel->FetchData(&spcdata,&M,&N,vector_enum);
    3939
    40         /*Partition edges if we are using P2 finite elements*/
    41         if(finite_element==P2Enum){
    42                 EdgesPartitioning(&my_edges,iomodel);
    43         }
    44 
     40        switch(finite_element){
     41                case P1Enum:
     42                        /*Nothing else to do*/
     43                        break;
     44                case P2Enum:
     45                        EdgesPartitioning(&my_edges,iomodel);
     46                        break;
     47                case P2xP1Enum:
     48                        EdgesPartitioning(&my_edges,iomodel);
     49                        break;
     50                default:
     51                        _error_("Finite element "<<EnumToStringx(finite_element)<<" not supported yet");
     52        }
     53
     54        count=0;
    4555        if(M==iomodel->numberofvertices){
    46                 /*Static constraint*/
    47                 count=0;
    48                 for (i=0;i<iomodel->numberofvertices;i++){
    49                         if((iomodel->my_vertices[i])){
    50                                 if (!xIsNan<IssmDouble>(spcdata[i])){
    51                                         constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,dof,spcdata[i],analysis_type));
    52                                         count++;
    53                                 }
    54                         }
    55                 }
    56                 if(finite_element==P2Enum){
    57                         for(i=0;i<iomodel->numberofedges;i++){
    58                                 if(my_edges[i]){
    59                                         v1 = iomodel->edges[3*i+0]-1;
    60                                         v2 = iomodel->edges[3*i+1]-1;
    61                                         if(!xIsNan<IssmDouble>(spcdata[v1]) && !xIsNan<IssmDouble>(spcdata[v2])){
    62                                                 constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+iomodel->numberofvertices+i+1,
    63                                                                                 dof,(spcdata[v1]+spcdata[v2])/2.,analysis_type));
    64                                                 count++;
    65                                         }
    66                                 }
    67                         }
     56                switch(finite_element){
     57                        case P1Enum:
     58                                for(i=0;i<iomodel->numberofvertices;i++){
     59                                        if((iomodel->my_vertices[i])){
     60                                                if (!xIsNan<IssmDouble>(spcdata[i])){
     61                                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,dof,spcdata[i],analysis_type));
     62                                                        count++;
     63                                                }
     64                                        }
     65                                }
     66                                break;
     67                        case P2Enum:
     68                                for(i=0;i<iomodel->numberofvertices;i++){
     69                                        if((iomodel->my_vertices[i])){
     70                                                if (!xIsNan<IssmDouble>(spcdata[i])){
     71                                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,dof,spcdata[i],analysis_type));
     72                                                        count++;
     73                                                }
     74                                        }
     75                                }
     76                                for(i=0;i<iomodel->numberofedges;i++){
     77                                        if(my_edges[i]){
     78                                                v1 = iomodel->edges[3*i+0]-1;
     79                                                v2 = iomodel->edges[3*i+1]-1;
     80                                                if(!xIsNan<IssmDouble>(spcdata[v1]) && !xIsNan<IssmDouble>(spcdata[v2])){
     81                                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+iomodel->numberofvertices+i+1,
     82                                                                                        dof,(spcdata[v1]+spcdata[v2])/2.,analysis_type));
     83                                                        count++;
     84                                                }
     85                                        }
     86                                }
     87                                break;
     88                        case P2xP1Enum:
     89                                for(i=0;i<iomodel->numberofvertices;i++){
     90                                        if((iomodel->my_vertices[i])){
     91                                                if (!xIsNan<IssmDouble>(spcdata[i])){
     92                                                        constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,dof,spcdata[i],analysis_type));
     93                                                        count++;
     94                                                }
     95                                        }
     96                                }
     97                                for(i=0;i<iomodel->numberofedges;i++){
     98                                        if(iomodel->edges[i*3+2]!=2){
     99                                                if(my_edges[i]){
     100                                                        v1 = iomodel->edges[3*i+0]-1;
     101                                                        v2 = iomodel->edges[3*i+1]-1;
     102                                                        if(!xIsNan<IssmDouble>(spcdata[v1]) && !xIsNan<IssmDouble>(spcdata[v2])){
     103                                                                constraints->AddObject(new SpcStatic(iomodel->constraintcounter+count+1,iomodel->nodecounter+iomodel->numberofvertices+i+1,
     104                                                                                                dof,(spcdata[v1]+spcdata[v2])/2.,analysis_type));
     105                                                                count++;
     106                                                        }
     107                                                }
     108                                        }
     109                                }
     110                                break;
     111                        default:
     112                                _error_("Finite element "<<EnumToStringx(finite_element)<<" not supported yet");
    68113                }
    69114        }
    … …  
    71116                /*transient: create transient SpcTransient objects. Same logic, except we need to retrieve
    72117                 * various times and values to initialize an SpcTransient object: */
    73                 count=0;
    74118
    75119                /*figure out times: */
    … …  
    77121                for(j=0;j<N;j++) times[j]=spcdata[(M-1)*N+j]*yts;
    78122
    79                 /*Create constraints: */
    80                 for (i=0;i<iomodel->numberofvertices;i++){
    81                         if((iomodel->my_vertices[i])){
    82 
    83                                 /*figure out times and values: */
    84                                 values=xNew<IssmDouble>(N);
    85                                 spcpresent=false;
    86                                 for(j=0;j<N;j++){
    87                                         values[j]=spcdata[i*N+j];
    88                                         if(!xIsNan<IssmDouble>(values[j]))spcpresent=true; //NaN means no spc by default
    89                                 }
    90 
    91                                 if(spcpresent){
    92                                         constraints->AddObject(new SpcTransient(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,dof,N,times,values,analysis_type));
    93                                         count++;
    94                                 }
    95                                 xDelete<IssmDouble>(values);
    96                         }
    97                 }
    98                 if(finite_element==P2Enum){
    99                         for(i=0;i<iomodel->numberofedges;i++){
    100                                 if(my_edges[i]){
    101                                         v1 = iomodel->edges[3*i+0]-1;
    102                                         v2 = iomodel->edges[3*i+1]-1;
    103                                         values=xNew<IssmDouble>(N);
    104                                         spcpresent=false;
    105                                         for(j=0;j<N;j++){
    106                                                 values[j]=(spcdata[v1*N+j]+spcdata[v2*N+j])/2.;
    107                                                 if(!xIsNan<IssmDouble>(values[j])) spcpresent=true; //NaN means no spc by default
    108                                         }
    109                                         if(spcpresent){
    110                                                 constraints->AddObject(new SpcTransient(iomodel->constraintcounter+count+1,iomodel->nodecounter+iomodel->numberofvertices+i+1,dof,
    111                                                                                 N,times,values,analysis_type));
    112                                                 count++;
    113                                         }
    114                                         xDelete<IssmDouble>(values);
    115                                 }
    116                         }
     123                switch(finite_element){
     124                        case P1Enum:
     125                                for(i=0;i<iomodel->numberofvertices;i++){
     126                                        if((iomodel->my_vertices[i])){
     127
     128                                                /*figure out times and values: */
     129                                                values=xNew<IssmDouble>(N);
     130                                                spcpresent=false;
     131                                                for(j=0;j<N;j++){
     132                                                        values[j]=spcdata[i*N+j];
     133                                                        if(!xIsNan<IssmDouble>(values[j]))spcpresent=true; //NaN means no spc by default
     134                                                }
     135
     136                                                if(spcpresent){
     137                                                        constraints->AddObject(new SpcTransient(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,dof,N,times,values,analysis_type));
     138                                                        count++;
     139                                                }
     140                                                xDelete<IssmDouble>(values);
     141                                        }
     142                                }
     143                                break;
     144                        case P2Enum:
     145                                for(i=0;i<iomodel->numberofvertices;i++){
     146                                        if((iomodel->my_vertices[i])){
     147
     148                                                /*figure out times and values: */
     149                                                values=xNew<IssmDouble>(N);
     150                                                spcpresent=false;
     151                                                for(j=0;j<N;j++){
     152                                                        values[j]=spcdata[i*N+j];
     153                                                        if(!xIsNan<IssmDouble>(values[j]))spcpresent=true; //NaN means no spc by default
     154                                                }
     155
     156                                                if(spcpresent){
     157                                                        constraints->AddObject(new SpcTransient(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,dof,N,times,values,analysis_type));
     158                                                        count++;
     159                                                }
     160                                                xDelete<IssmDouble>(values);
     161                                        }
     162                                }
     163                                for(i=0;i<iomodel->numberofedges;i++){
     164                                        if(my_edges[i]){
     165                                                v1 = iomodel->edges[3*i+0]-1;
     166                                                v2 = iomodel->edges[3*i+1]-1;
     167                                                values=xNew<IssmDouble>(N);
     168                                                spcpresent=false;
     169                                                for(j=0;j<N;j++){
     170                                                        values[j]=(spcdata[v1*N+j]+spcdata[v2*N+j])/2.;
     171                                                        if(!xIsNan<IssmDouble>(values[j])) spcpresent=true; //NaN means no spc by default
     172                                                }
     173                                                if(spcpresent){
     174                                                        constraints->AddObject(new SpcTransient(iomodel->constraintcounter+count+1,iomodel->nodecounter+iomodel->numberofvertices+i+1,dof,
     175                                                                                        N,times,values,analysis_type));
     176                                                        count++;
     177                                                }
     178                                                xDelete<IssmDouble>(values);
     179                                        }
     180                                }
     181                                break;
     182                        case P2xP1Enum:
     183                                for(i=0;i<iomodel->numberofvertices;i++){
     184                                        if((iomodel->my_vertices[i])){
     185
     186                                                /*figure out times and values: */
     187                                                values=xNew<IssmDouble>(N);
     188                                                spcpresent=false;
     189                                                for(j=0;j<N;j++){
     190                                                        values[j]=spcdata[i*N+j];
     191                                                        if(!xIsNan<IssmDouble>(values[j]))spcpresent=true; //NaN means no spc by default
     192                                                }
     193
     194                                                if(spcpresent){
     195                                                        constraints->AddObject(new SpcTransient(iomodel->constraintcounter+count+1,iomodel->nodecounter+i+1,dof,N,times,values,analysis_type));
     196                                                        count++;
     197                                                }
     198                                                xDelete<IssmDouble>(values);
     199                                        }
     200                                }
     201                                for(i=0;i<iomodel->numberofedges;i++){
     202                                        if(iomodel->edges[i*3+2]!=2){
     203                                                if(my_edges[i]){
     204                                                        v1 = iomodel->edges[3*i+0]-1;
     205                                                        v2 = iomodel->edges[3*i+1]-1;
     206                                                        values=xNew<IssmDouble>(N);
     207                                                        spcpresent=false;
     208                                                        for(j=0;j<N;j++){
     209                                                                values[j]=(spcdata[v1*N+j]+spcdata[v2*N+j])/2.;
     210                                                                if(!xIsNan<IssmDouble>(values[j])) spcpresent=true; //NaN means no spc by default
     211                                                        }
     212                                                        if(spcpresent){
     213                                                                constraints->AddObject(new SpcTransient(iomodel->constraintcounter+count+1,iomodel->nodecounter+iomodel->numberofvertices+i+1,dof,
     214                                                                                                N,times,values,analysis_type));
     215                                                                count++;
     216                                                        }
     217                                                        xDelete<IssmDouble>(values);
     218                                                }
     219                                        }
     220                                }
     221                                break;
     222                        default:
     223                                _error_("Finite element "<<EnumToStringx(finite_element)<<" not supported yet");
    117224                }
    118225        }
Note: See TracChangeset for help on using the changeset viewer.