Changeset 4013


Ignore:
Timestamp:
Jun 3, 2010, 12:09:20 PM (16 years ago)
Author:
Eric.Larour
Message:

New parallel diagnostic solution

Location:
issm/trunk/src/c
Files:
1 added
14 edited

Legend:

Unmodified
Added
Removed
  • issm/trunk/src/c/EnumDefinitions/EnumDefinitions.h

    r3982 r4013  
    4646        //slope
    4747        SurfaceAnalysisEnum,
    48         SlopecomputeAnalysisEnum,
     48        SlopeComputeAnalysisEnum,
    4949        BedXAnalysisEnum,
    5050        BedYAnalysisEnum,
    … …  
    154154        ArtDiffEnum,
    155155        BedEnum,
    156         BedSlopexEnum,
    157         BedSlopeyEnum,
     156        BedSlopeXEnum,
     157        BedSlopeYEnum,
    158158        BoundaryEnum,
    159159        CmMaxDmpSlopeEnum,
    … …  
    204204        SurfaceAreaEnum,
    205205        SurfaceEnum,
    206         SurfaceSlopexEnum,
    207         SurfaceSlopeyEnum,
     206        SurfaceSlopeXEnum,
     207        SurfaceSlopeYEnum,
    208208        TemperatureEnum,
    209209        TemperatureAverageEnum,
  • issm/trunk/src/c/Makefile.am

    r4002 r4013  
    951951                                        ./solutions/transient_core_3d.cpp\
    952952                                        ./solutions/steadystate_core.cpp\
     953                                        ./solutions/ResetBoundaryConditions.cpp\
     954                                        ./solutions/solve_linear.cpp\
    953955                                        ./modules/Bamgx/Bamgx.cpp\
    954956                                        ./modules/Bamgx/Bamgx.h\
  • issm/trunk/src/c/modules/GetSolutionFromInputsx/GetSolutionFromInputsx.h

    r3913 r4013  
    1010
    1111/* local prototypes: */
    12 void    GetSolutionFromInputsx( Vec* psolution, DataSet* elements,DataSet* nodes, DataSet* vertices,DataSet* loads, DataSet* materials,  Parameters* parameters, int analysis_type, int sub_analysis_type);
     12void GetSolutionFromInputsx( Vec* psolution, DataSet* elements,DataSet* nodes, DataSet* vertices,DataSet* loads, DataSet* materials,  Parameters* parameters, int analysis_type, int sub_analysis_type);
    1313
    1414#endif  /* _GETSOLUTIONFROMINPUTSXX_H */
  • issm/trunk/src/c/modules/OutputResultsx/OutputResultsx.cpp

    r3938 r4013  
    3535                pfclose(fid,filename);
    3636        }
     37
     38
     39        /*output if we have rifts: */
     40        if(numrifts){
     41                OutputRiftsx( &riftproperties,fem_dh->loads,numrifts);
     42                results->AddObject(new Result(results->Size()+1,0,1,"riftproperties",riftproperties));
     43        }
     44
     45
    3746}
  • issm/trunk/src/c/objects/Elements/Beam.cpp

    r4011 r4013  
    119119/*FUNCTION Beam::IsInput{{{1*/
    120120bool Beam::IsInput(int name){
    121         if (name==SurfaceSlopexEnum ||
    122                                 name==SurfaceSlopeyEnum){
     121        if (name==SurfaceSlopeXEnum ||
     122                                name==SurfaceSlopeYEnum){
    123123                return true;
    124124        }
    … …  
    412412
    413413        //recover extra inputs
    414         inputs->GetParameterValue(&slope[0],&gauss1[0],SurfaceSlopexEnum);
    415         inputs->GetParameterValue(&slope[1],&gauss1[0],SurfaceSlopeyEnum);
     414        inputs->GetParameterValue(&slope[0],&gauss1[0],SurfaceSlopeXEnum);
     415        inputs->GetParameterValue(&slope[1],&gauss1[0],SurfaceSlopeYEnum);
    416416        inputs->GetParameterValue(&surface,&gauss1[0],SurfaceEnum);
    417417        inputs->GetParameterValue(&thickness,&gauss1[0],ThicknessEnum);
  • issm/trunk/src/c/objects/Elements/Penta.cpp

    r4011 r4013  
    577577
    578578        }
    579         else if (analysis_type==SlopecomputeAnalysisEnum){
     579        else if (analysis_type==SlopeComputeAnalysisEnum){
    580580
    581581                UpdateInputsFromSolutionSlopeCompute( solution,analysis_type,sub_analysis_type);
    … …  
    813813                                name==SurfaceEnum ||
    814814                                name==BedEnum ||
    815                                 name==SurfaceSlopexEnum ||
    816                                 name==SurfaceSlopeyEnum ||
     815                                name==SurfaceSlopeXEnum ||
     816                                name==SurfaceSlopeYEnum ||
    817817                                name==MeltingRateEnum ||
    818818                                name==AccumulationRateEnum ||
  • issm/trunk/src/c/objects/Elements/Sing.cpp

    r4011 r4013  
    117117/*FUNCTION Sing::IsInput{{{1*/
    118118bool Sing::IsInput(int name){
    119         if (name==SurfaceSlopexEnum ||
    120                                 name==SurfaceSlopeyEnum){
     119        if (name==SurfaceSlopeXEnum ||
     120                                name==SurfaceSlopeYEnum){
    121121                return true;
    122122        }
    … …  
    324324        matice=(Matice*)hmatice.delivers();
    325325
    326         inputs->GetParameterValue(&slope[0],SurfaceSlopexEnum);
    327         inputs->GetParameterValue(&slope[1],SurfaceSlopeyEnum);
     326        inputs->GetParameterValue(&slope[0],SurfaceSlopeXEnum);
     327        inputs->GetParameterValue(&slope[1],SurfaceSlopeYEnum);
    328328
    329329        GetDofList(&doflist[0],&numberofdofspernode);
  • issm/trunk/src/c/objects/Elements/Tria.cpp

    r4011 r4013  
    501501                else ISSMERROR("%s%i%s\n","sub_analysis: ",sub_analysis_type," not supported yet");
    502502        }
    503         else if (analysis_type==SlopecomputeAnalysisEnum){
     503        else if (analysis_type==SlopeComputeAnalysisEnum){
    504504                UpdateInputsFromSolutionSlopeCompute( solution,analysis_type,sub_analysis_type);
    505505        }
    … …  
    693693/*FUNCTION Tria::IsInput{{{1*/
    694694bool Tria::IsInput(int name){
    695         if (name==SurfaceSlopexEnum ||
    696                                 name==SurfaceSlopeyEnum){
     695        if (name==SurfaceSlopeXEnum ||
     696                                name==SurfaceSlopeYEnum){
    697697                return true;
    698698        }
    … …  
    878878                else ISSMERROR("%s%i%s\n","sub_analysis: ",sub_analysis_type," not supported yet");
    879879        }
    880         else if (analysis_type==SlopecomputeAnalysisEnum){
     880        else if (analysis_type==SlopeComputeAnalysisEnum){
    881881                CreateKMatrixSlopeCompute( Kgg,analysis_type,sub_analysis_type);
    882882        }
  • issm/trunk/src/c/solutions/diagnostic.cpp

    r4004 r4013  
    3131
    3232        /*Results: */
    33         Results* results=NULL;
    34        
    3533        bool waitonlock=false;
    3634       
    … …  
    8381                        _printf_("call computational core:\n");
    8482                        MPI_Barrier(MPI_COMM_WORLD); start_core=MPI_Wtime( );
    85                         results=diagnostic_core(femmodel);
     83                        iagnostic_core(femmodel);
    8684                        MPI_Barrier(MPI_COMM_WORLD); finish_core=MPI_Wtime( );
    8785
    … …  
    9189                        _printf_("call computational core:\n");
    9290                        MPI_Barrier(MPI_COMM_WORLD); start_core=MPI_Wtime( );
    93                         results=control_core(femmodel);
     91                        ontrol_core(femmodel);
    9492                        MPI_Barrier(MPI_COMM_WORLD); finish_core=MPI_Wtime( );
    9593
    … …  
    9795
    9896                _printf_("write results to disk:\n");
    99                 OutputResults(results,outputfilename);
     97                OutputResults(outputfilename);
    10098        }
    10199        else{
    … …  
    119117        /*Free ressources */
    120118        delete femmodel;
    121         delete results;
    122119
    123120        /*Get finish time and close*/
  • issm/trunk/src/c/solutions/diagnostic_core.cpp

    r4012 r4013  
    1111#include "../include/include.h"
    1212
    13 Results* diagnostic_core(FemModel* femmodel){
    14 
    15         extern int my_rank;
    16         int        dummy;
    17 
    18         /*output: */
    19         Results* results=NULL;
    20         Result*  result=NULL;
     13void* diagnostic_core(FemModel* femmodel){
    2114
    2215        /*solutions: */
    2316        Vec ug=NULL;
    24         Vec ug_horiz=NULL;
    25         Vec ug_vert=NULL;
    26         Vec ug_stokes=NULL;
    27         Vec pg=NULL;
    28         Vec riftproperties=NULL;
    29         double* u_g_initial=NULL;
    30         double* vx=NULL;
    31         double* vy=NULL;
    32         double* vz=NULL;
    3317
    34         /*flags: */
    35         int verbose=0;
     18        /*parameters: */
     19        int  verbose=0;
    3620        bool qmu_analysis=false;
    37         int dim=-1;
     21        int  dim=-1;
    3822        bool ishutter=false;
    3923        bool ismacayealpattyn=false;
    4024        bool isstokes=false;
    41         int numberofdofspernode_sl;
    42         int numberofdofspernode_dh;
    43         int numberofdofspernode_ds;
    44         int numberofnodes;
    45         int numrifts=0;
    46 
    47         /*slopes: */
    48         Vec surfaceslopex=NULL;
    49         Vec surfaceslopey=NULL;
    50         Vec bedslopex=NULL;
    51         Vec bedslopey=NULL;
    52 
    5325        double stokesreconditioning;
    54 
    55         /*dof recovery: */
    56         int dof01[2]={0,1};
    57         int dof2[1]={2};
    58         int dof012[3]={0,1,2};
    59         int dof3[1]={3};
    60         double* dofset=NULL;
    61 
    62         //initialize results
    63         results=new Results();
    6426
    6527        //first recover parameters needed to drive the solution
    … …  
    6830        femmodel->parameters->FindParam(&ishutter,IsHutterEnum);
    6931        femmodel->parameters->FindParam(&ismacayealpattyn,IsMacAyealPattynEnum);
    70         femmodel->parameters->FindParam(&numberofnodes,NumberOfNodesEnum);
    7132        femmodel->parameters->FindParam(&isstokes,IsStokesEnum);
    7233        femmodel->parameters->FindParam(&stokesreconditioning,StokesReconditioningEnum);
    73         femmodel->parameters->FindParam(&numrifts,NumRiftsEnum);
    7434        femmodel->parameters->FindParam(&qmu_analysis,QmuAnalysisEnum);
    75 
    76         //specific parameters for specific models
    77         femmodel->parameters->FindParam(&numberofdofspernode_dh,NumberOfDofsPerNodeEnum);
    78         femmodel->parameters->FindParam(&numberofdofspernode_sl,NumberOfDofsPerNodeEnum);
    79         femmodel->parameters->FindParam(&numberofdofspernode_ds,NumberOfDofsPerNodeEnum);
    8035
    8136        /*for qmu analysis, reinitialize velocity so that fake sensitivities do not show up as a result of a different restart of the convergence at each trial.*/
    … …  
    8742
    8843        /*Compute slopes: */
    89         slope_core(femmodel,SurfaceXAnalysisEnum);
    90         slope_core(femmodel,SurfaceYAnalysisEnum);
    91         slope_core(femmodel,BedSlopeXAnalysisEnum);
    92         slope_core(femmodel,BedSlopeYAnalysisEnum);
     44        if(ishutter){
     45                slope_core(femmodel,SurfaceXAnalysisEnum);
     46                slope_core(femmodel,SurfaceYAnalysisEnum);
     47        }
     48        if(isstokes){
     49                slope_core(femmodel,BedSlopeXAnalysisEnum);
     50                slope_core(femmodel,BedSlopeYAnalysisEnum);
     51        }
    9352               
    9453       
    … …  
    9655                       
    9756                if(verbose)_printf_("%s\n"," computing hutter velocities...");
    98                 diagnostic_core_linear(&ug,fem_dhu,DiagnosticAnalysisEnum,HutterAnalysisEnum);
    99 
    100                 if(verbose)_printf_("%s\n"," computing pressure according to MacAyeal...");
    101                 ComputePressurex(&pg,fem_dhu->elements,fem_dhu->nodes, fem_dhu->vertices,fem_dhu->loads,fem_dhu->materials,fem_dhu->parameters,DiagnosticAnalysisEnum,HutterAnalysisEnum);
    102 
    103                 if(verbose)_printf_("%s\n"," update boundary conditions for macyeal pattyn using hutter results...");
    104                 if (ismacayealpattyn){
    105                         VecFree(&fem_dh->yg->vector); VecFree(&fem_dh->ys);
    106                         VecDuplicatePatch(&fem_dh->yg->vector,ug);
    107                         Reducevectorgtosx(&fem_dh->ys,fem_dh->yg->vector,fem_dh->nodesets);
    108                 }
    109 
     57                solve_linear(&ug,femmodel,DiagnosticAnalysisEnum,HutterAnalysisEnum);
     58               
     59                if (ismacayealpattyn) ResetBoundaryConditions(femmodel,DiagnosticHorizAnalysisEnum,ug);
     60                VecFree(&ug);
    11061        }
    11162
    … …  
    11364               
    11465                if(verbose)_printf_("%s\n"," computing horizontal velocities...");
    115                 diagnostic_core_nonlinear(&ug,NULL,NULL,fem_dh->loads,fem_dh,DiagnosticAnalysisEnum,HorizAnalysisEnum);
    116 
    117                 if(dim==2){
    118                         if(verbose)_printf_("%s\n"," computing pressure according to MacAyeal...");
    119                         ComputePressurex(&pg,fem_dh->elements,fem_dh->nodes, fem_dh->vertices,fem_dh->loads,fem_dh->materials,fem_dh->parameters,DiagnosticAnalysisEnum,HorizAnalysisEnum);
    120                 }
    121 
     66                diagnostic_core_nonlinear(NULL,NULL,NULL,femmodel,false,DiagnosticAnalysisEnum,HorizAnalysisEnum); //false means we expose loads to changes inside the solution
    12267        }
    12368       
    … …  
    12570        if (dim==3){
    12671
    127                 if(verbose)_printf_("%s\n"," extruding horizontal velocities...");
    128                 FieldExtrudex(ug,fem_dh->elements,fem_dh->nodes, fem_dh->vertices,fem_dh->loads,fem_dh->materials,fem_dh->parameters,"velocity",1);
    129                 SplitSolutionVectorx(ug,numberofnodes,numberofdofspernode_dh,&vx,&vy);
    130                 model->UpdateInputsFromVector(vx,VxEnum,VertexEnum);
    131                 model->UpdateInputsFromVector(vy,VyEnum,VertexEnum);
    132                
    13372                if(verbose)_printf_("%s\n"," computing vertical velocities...");
    134                 diagnostic_core_linear(&ug_vert,fem_dv,DiagnosticAnalysisEnum,VertAnalysisEnum);
    135                 model->UpdateInputsFromVector(ug_vert,VzEnum,VertexEnum);
     73                solve_linear(NULL,femmodel,DiagnosticAnalysisEnum,VertAnalysisEnum);
    13674
    137                 if(verbose)_printf_("%s\n"," computing pressure according to Pattyn...");
    138                 ComputePressurex(&pg,fem_dh->elements, fem_dh->nodes, fem_dh->vertices,fem_dh->loads,  fem_dh->materials,fem_dh->parameters,DiagnosticAnalysisEnum,HorizAnalysisEnum);
    139                 model->UpdateInputsFromVector(pg,PressureEnum,VertexEnum);
    140                
    14175                if (isstokes){
    14276
    143                         //"recondition" pressure
    144                         VecScale(pg,1.0/stokesreconditioning);
     77                        //"recondition" pressure computed previously:
     78                        DuplicateInputx(femmodel,PressureEnum,PressureStokesEnum); ScaleInputx(femmmodel,PressureStokesEnum,1.0/stokesreconditioning);
    14579
    14680                        if(verbose)_printf_("%s\n"," update boundary conditions for stokes using velocities previously computed...");
    147                         xfree((void**)&dofset);dofset=dofsetgen(2,dof01,4,numberofnodes*4); VecMerge(fem_ds->yg->vector,ug,dofset,2*numberofnodes);
    148                         xfree((void**)&dofset);dofset=dofsetgen(1,dof2,4,numberofnodes*4); VecMerge(fem_ds->yg->vector,ug_vert,dofset,1*numberofnodes);
    149                         VecFree(&fem_ds->ys);
    150                         Reducevectorgtosx(&fem_ds->ys,fem_ds->yg->vector,fem_ds->nodesets);
     81                        GetSolutionFromInputsx( &ug, femmodel->elements,femmodel->nodes, femmodel->vertices,femmodel->loads, femmodel->materials,  femmodel->parameters, DiagnosticAnalysisEnum,StokesAnalysisEnum);
     82                        ResetBoundaryConditions(femmodel,DiagnosticStokesAnalysisEnum,ug);
     83                        VecFree(&ug);
    15184
    15285                        if(verbose)_printf_("%s\n"," computing stokes velocities and pressure ...");
    153                         VecFree(&ug);
    154                         diagnostic_core_nonlinear(&ug,NULL,NULL,NULL,fem_ds,DiagnosticAnalysisEnum,StokesAnalysisEnum);
    155                
    156                         //decondition" pressure
    157                         VecFree(&pg);   
    158                         xfree((void**)&dofset);dofset=dofsetgen(1,dof3,4,numberofnodes*4); VecPartition(&pg, ug, dofset, numberofnodes*1);
    159                         VecScale(pg,stokesreconditioning);
     86                        diagnostic_core_nonlinear(NULL,NULL,NULL,NULL,fem_ds,DiagnosticAnalysisEnum,StokesAnalysisEnum);
    16087                }
    16188        }
    16289       
    163         /*Plug results into output dataset: */
    164         if(dim==2){
    165                 if(ismacayealpattyn){
    166                         InputToResultx(&result,fem_dh->elements,fem_dh->nodes,fem_dh->vertices, fem_dh->loads, fem_dh->materials,fem_dh->parameters,VxEnum,results->Size()+1,0,1); results->AddObject(result);
    167                         InputToResultx(&result,fem_dh->elements,fem_dh->nodes,fem_dh->vertices, fem_dh->loads, fem_dh->materials,fem_dh->parameters,VyEnum,results->Size()+1,0,1); results->AddObject(result);
    168                         InputToResultx(&result,fem_dh->elements,fem_dh->nodes,fem_dh->vertices, fem_dh->loads, fem_dh->materials,fem_dh->parameters,PressureEnum,results->Size()+1,0,1); results->AddObject(result);
    169                 }
    170                 else{
    171                         InputToResultx(&result,fem_dhu->elements,fem_dhu->nodes,fem_dhu->vertices, fem_dhu->loads, fem_dhu->materials,fem_dhu->parameters,VxEnum,results->Size()+1,0,1); results->AddObject(result);
    172                         InputToResultx(&result,fem_dhu->elements,fem_dhu->nodes,fem_dhu->vertices, fem_dhu->loads, fem_dhu->materials,fem_dhu->parameters,VyEnum,results->Size()+1,0,1); results->AddObject(result);
    173                         InputToResultx(&result,fem_dhu->elements,fem_dhu->nodes,fem_dhu->vertices, fem_dhu->loads, fem_dhu->materials,fem_dhu->parameters,PressureEnum,results->Size()+1,0,1); results->AddObject(result);
    174                 }
    175         }
    176         else{
    177                 if(isstokes){
    178                         InputToResultx(&result,fem_ds->elements,fem_ds->nodes,fem_ds->vertices, fem_ds->loads, fem_ds->materials,fem_ds->parameters,VxEnum,results->Size()+1,0,1); results->AddObject(result);
    179                         InputToResultx(&result,fem_ds->elements,fem_ds->nodes,fem_ds->vertices, fem_ds->loads, fem_ds->materials,fem_ds->parameters,VyEnum,results->Size()+1,0,1); results->AddObject(result);
    180                         InputToResultx(&result,fem_ds->elements,fem_ds->nodes,fem_ds->vertices, fem_ds->loads, fem_ds->materials,fem_ds->parameters,VzEnum,results->Size()+1,0,1); results->AddObject(result);
    181                         InputToResultx(&result,fem_ds->elements,fem_ds->nodes,fem_ds->vertices, fem_ds->loads, fem_ds->materials,fem_ds->parameters,PressureEnum,results->Size()+1,0,1); results->AddObject(result);
    182                 }
    183                 else{
    184                         if(ismacayealpattyn){
    185                                 InputToResultx(&result,fem_dh->elements,fem_dh->nodes,fem_dh->vertices, fem_dh->loads, fem_dh->materials,fem_dh->parameters,VxEnum,results->Size()+1,0,1); results->AddObject(result);
    186                                 InputToResultx(&result,fem_dh->elements,fem_dh->nodes,fem_dh->vertices, fem_dh->loads, fem_dh->materials,fem_dh->parameters,VyEnum,results->Size()+1,0,1); results->AddObject(result);
    187                                 InputToResultx(&result,fem_dh->elements,fem_dh->nodes,fem_dh->vertices, fem_dh->loads, fem_dh->materials,fem_dh->parameters,VzEnum,results->Size()+1,0,1); results->AddObject(result);
    188                                 InputToResultx(&result,fem_dh->elements,fem_dh->nodes,fem_dh->vertices, fem_dh->loads, fem_dh->materials,fem_dh->parameters,PressureEnum,results->Size()+1,0,1); results->AddObject(result);
    189                         }
    190                         else{
    191                                 InputToResultx(&result,fem_dhu->elements,fem_dhu->nodes,fem_dhu->vertices, fem_dhu->loads, fem_dhu->materials,fem_dhu->parameters,VxEnum,results->Size()+1,0,1); results->AddObject(result);
    192                                 InputToResultx(&result,fem_dhu->elements,fem_dhu->nodes,fem_dhu->vertices, fem_dhu->loads, fem_dhu->materials,fem_dhu->parameters,VyEnum,results->Size()+1,0,1); results->AddObject(result);
    193                                 InputToResultx(&result,fem_dhu->elements,fem_dhu->nodes,fem_dhu->vertices, fem_dhu->loads, fem_dhu->materials,fem_dhu->parameters,VzEnum,results->Size()+1,0,1); results->AddObject(result);
    194                                 InputToResultx(&result,fem_dhu->elements,fem_dhu->nodes,fem_dhu->vertices, fem_dhu->loads, fem_dhu->materials,fem_dhu->parameters,PressureEnum,results->Size()+1,0,1); results->AddObject(result);
    195                         }
    196 
    197                 }
    198         }
    199        
    200         results->AddObject(new StringResult(results->Size()+1,AnalysisTypeEnum,0,1,EnumAsString(DiagnosticAnalysisEnum)));
    201 
    202         /*output if we have rifts: */
    203         if(numrifts){
    204                 OutputRiftsx( &riftproperties,fem_dh->loads,numrifts);
    205                 results->AddObject(new Result(results->Size()+1,0,1,"riftproperties",riftproperties));
    206         }
    207 
    20890        /*Free ressources: */
    20991        VecFree(&ug);
    210         VecFree(&pg);
    211         xfree((void**)&dofset);
    212         xfree((void**)&u_g_initial);
    213        
    214         VecFree(&surfaceslopex);
    215         VecFree(&surfaceslopey);
    216         VecFree(&bedslopex);
    217         VecFree(&bedslopey);
    218 
    219         /*return: */
    220         return results;
    22192}
  • issm/trunk/src/c/solutions/diagnostic_core_nonlinear.cpp

    r3922 r4013  
    99#include "./solutions.h"
    1010
    11 void diagnostic_core_nonlinear(Vec* pug,Mat* pKff0,Mat* pKfs0, DataSet* input_loads,FemModel* fem,int analysis_type,int sub_analysis_type){
     11void diagnostic_core_nonlinear(Vec* pug,Mat* pKff0,Mat* pKfs0, FemModel* fem,bool conserve_loads, int analysis_type,int sub_analysis_type){
    1212
    1313
    … …  
    3838        int verbose=0;
    3939
     40        //Set active analysis_type:
     41        fem->SetActiveAnalysis(analysis_type);
     42
    4043        /*Recover parameters: */
    4144        kflag=1; pflag=1;
    … …  
    4750       
    4851        /*Were loads requested as output? : */
    49         if(!input_loads){
    50                 loads=fem->loads->Copy(); //we don't want to clobber loads, as they are not needed in output.
     52        if(conserve_loads){
     53                loads=fem->loads->Copy(); //protect loads from being modified by the solution
    5154        }
    5255        else{
    53                 loads=input_loads; //we are going to modify the loads!
     56                loads=fem->loads; //modify loads  in this solution
    5457        }
    5558
    … …  
    130133                }
    131134        }
     135       
     136        /*extrude if we are in 3D: */
     137        if (dim==3){
     138                if(verbose)_printf_("%s\n","extruding velocity and pressure in 3d...");
     139                InputsExtrudex( fem->elements,fem->nodes,fem->vertices,fem->loads,fem->materials,fem->parameters,VxEnum,0);
     140                InputsExtrudex( fem->elements,fem->nodes,fem->vertices,fem->loads,fem->materials,fem->parameters,VyEnum,0);
     141        }
    132142
    133143        //more output might be needed, when running in control
    … …  
    143153
    144154        /*Delete loads only if no ouput was requested: */
    145         if(!input_loads)delete loads;
     155        if(!conserve_loads)delete loads;
    146156
    147157        /*clean up*/
    … …  
    152162       
    153163        /*Assign output pointers: */
    154         *pug=ug;
     164        if(pug)*pug=ug;
     165        else VecFree(&ug);
     166
    155167        if(pKff0)*pKff0=Kff;
    156168        if(pKfs0)*pKfs0=Kfs;
  • issm/trunk/src/c/solutions/slope_core.cpp

    r4012 r4013  
    2828        if(verbose)_printf_("%s\n","computing slope...");
    2929
    30         /*Early return possible? */
    31         if(!ishutter && AnalysisEnum==SurfaceAnalysisEnum) return;
    32         if(!isstokes && AnalysisEnum==BedAnalysisEnum) return;
    33 
    3430        /*Call on core computations: */
    3531        solver_linear(&slope,femmodel,SlopecomputeAnalysisEnum,sub_analysis_type);
  • issm/trunk/src/c/solutions/solutions.h

    r4012 r4013  
    1414/*cores: */
    1515Results* gradjcompute_core(Model* model);
    16 Results* diagnostic_core(Model* model);
     16void* diagnostic_core(Model* model);
    1717Results* prognostic_core(Model* model);
    1818Results* prognostic2_core(Model* model);
    … …  
    2727Results* transient_core_3d(Model* model);
    2828Results* thermal_core(Model* model);
     29void slope_core(FemModel* fem,int sub_analysis_type);
    2930
    3031/*computational cores: */
    3132void thermal_core_nonlinear(Vec* ptg,double* pmelting_offset,FemModel* fem,int analysis_type,int sub_analysis_type);
    32 void diagnostic_core_nonlinear(Vec* pug,Mat* pK_ff0,Mat* pK_fs0, DataSet* loads, FemModel* fem,int analysis_type,int sub_analysis_type);
     33void diagnostic_core_nonlinear(Vec* pug,Mat* pK_ff0,Mat* pK_fs0, FemModel* fem,bool conserve_loads,int analysis_type,int sub_analysis_type);
    3334void diagnostic_core_linear(Vec* ppg,FemModel* fem,int  analysis_type,int sub_analysis_type);
    3435void solver_linear(Vec* ppg,FemModel* fem,int  analysis_type,int sub_analysis_type);
    35 void slope_core(FemModel* fem,int sub_analysis_type);
    3636
    3737
    … …  
    5656void CreateFemModel(FemModel* femmodel,ConstDataHandle MODEL,int analysis_type,int sub_analysis_type);
    5757//int BatchDebug(Mat* Kgg,Vec* pg,FemModel* femmodel,char* filename);
     58void ResetBoundaryConditions(FemModel* femmodel, int analysis_type, Vec ug);
    5859
    5960#endif
  • issm/trunk/src/c/solutions/solver_linear.cpp

    r4012 r4013  
    6161        xfree((void**)&solver_string);
    6262
    63         /*Assign output pointers:*/
    64         *pug=ug;
     63        /*Assign output pointers if requested:*/
     64        if (pug) *pug=ug;
    6565}
Note: See TracChangeset for help on using the changeset viewer.