source: issm/trunk-jpl/src/c/analyses/HydrologyShreveAnalysis.cpp@ 16761

#include "./HydrologyShreveAnalysis.h"
#include "../toolkits/toolkits.h"
#include "../classes/classes.h"
#include "../shared/shared.h"
#include "../modules/modules.h"

/*Model processing*/
int  HydrologyShreveAnalysis::DofsPerNode(int** doflist,int meshtype,int approximation){/*{{{*/
        return 1;
}/*}}}*/
void HydrologyShreveAnalysis::UpdateParameters(Parameters* parameters,IoModel* iomodel,int solution_enum,int analysis_enum){/*{{{*/

        /*retrieve some parameters: */
        int  hydrology_model;
        iomodel->Constant(&hydrology_model,HydrologyModelEnum);

        /*Now, do we really want Shreve?*/
        if(hydrology_model!=HydrologyshreveEnum) return;

        parameters->AddObject(new IntParam(HydrologyModelEnum,hydrology_model));
        parameters->AddObject(iomodel->CopyConstantObject(HydrologyshreveStabilizationEnum));

}/*}}}*/
void HydrologyShreveAnalysis::UpdateElements(Elements* elements,IoModel* iomodel,int analysis_counter,int analysis_type){/*{{{*/

        /*Fetch data needed: */
        int    hydrology_model;
        iomodel->Constant(&hydrology_model,HydrologyModelEnum);

        /*Now, do we really want Shreve?*/
        if(hydrology_model!=HydrologyshreveEnum) return;

        /*Update elements: */
        int counter=0;
        for(int i=0;i<iomodel->numberofelements;i++){
                if(iomodel->my_elements[i]){
                        Element* element=(Element*)elements->GetObjectByOffset(counter);
                        element->Update(i,iomodel,analysis_counter,analysis_type,P1Enum);
                        counter++;
                }
        }

        iomodel->FetchDataToInput(elements,ThicknessEnum);
        iomodel->FetchDataToInput(elements,SurfaceEnum);
        iomodel->FetchDataToInput(elements,BedEnum);
        iomodel->FetchDataToInput(elements,MeshElementonbedEnum);
        iomodel->FetchDataToInput(elements,MeshElementonsurfaceEnum);
        iomodel->FetchDataToInput(elements,MaskIceLevelsetEnum);
        iomodel->FetchDataToInput(elements,MaskGroundediceLevelsetEnum);
        iomodel->FetchDataToInput(elements,BasalforcingsMeltingRateEnum);
        iomodel->FetchDataToInput(elements,WatercolumnEnum);

        elements->InputDuplicate(WatercolumnEnum,WaterColumnOldEnum);
}/*}}}*/
void HydrologyShreveAnalysis::CreateNodes(Nodes* nodes,IoModel* iomodel){/*{{{*/

        /*Fetch parameters: */
        int  hydrology_model;
        iomodel->Constant(&hydrology_model,HydrologyModelEnum);

        /*Now, do we really want Shreve?*/
        if(hydrology_model!=HydrologyshreveEnum) return;

        if(iomodel->meshtype==Mesh3DEnum) iomodel->FetchData(2,MeshVertexonbedEnum,MeshVertexonsurfaceEnum);
        ::CreateNodes(nodes,iomodel,HydrologyShreveAnalysisEnum,P1Enum);
        iomodel->DeleteData(2,MeshVertexonbedEnum,MeshVertexonsurfaceEnum);
}/*}}}*/
void HydrologyShreveAnalysis::CreateConstraints(Constraints* constraints,IoModel* iomodel){/*{{{*/

        /*retrieve some parameters: */
        int          hydrology_model;
        iomodel->Constant(&hydrology_model,HydrologyModelEnum);

        if(hydrology_model!=HydrologyshreveEnum) return;

        IoModelToConstraintsx(constraints,iomodel,HydrologyshreveSpcwatercolumnEnum,HydrologyShreveAnalysisEnum,P1Enum);

}/*}}}*/
void HydrologyShreveAnalysis::CreateLoads(Loads* loads, IoModel* iomodel){/*{{{*/
        /*No loads*/
}/*}}}*/

/*Numerics*/
void HydrologyShreveAnalysis::GetSolutionFromInputs(Vector<IssmDouble>* solution,Element* element){/*{{{*/
        element->GetSolutionFromInputsOneDof(solution,WatercolumnEnum);
}/*}}}*/
void HydrologyShreveAnalysis::InputUpdateFromSolution(IssmDouble* solution,Element* element){/*{{{*/

        /*Intermediary*/
        int* doflist = NULL;

        /*Fetch number of nodes for this finite element*/
        int numnodes = element->GetNumberOfNodes();

        /*Fetch dof list and allocate solution vector*/
        element->GetDofList(&doflist,NoneApproximationEnum,GsetEnum);
        IssmDouble* values = xNew<IssmDouble>(numnodes);

        /*Use the dof list to index into the solution vector: */
        for(int i=0;i<numnodes;i++){
                values[i]=solution[doflist[i]];
                if(xIsNan<IssmDouble>(values[i])) _error_("NaN found in solution vector");
                if (values[i]<10e-10) values[i]=10e-10; //correcting the water column to positive values
        }

        /*Add input to the element: */
        element->AddInput(WatercolumnEnum,values,P1Enum);

        /*Free ressources:*/
        xDelete<IssmDouble>(values);
        xDelete<int>(doflist);
}/*}}}*/