Index: /issm/trunk-jpl/src/c/Makefile.am =================================================================== --- /issm/trunk-jpl/src/c/Makefile.am (revision 16381) +++ /issm/trunk-jpl/src/c/Makefile.am (revision 16382) @@ -92,4 +92,6 @@ ./classes/Inputs/Input.h\ ./classes/Inputs/InputLocal.h\ + ./classes/Inputs/SegInput.h\ + ./classes/Inputs/SegInput.cpp\ ./classes/Inputs/TriaInput.h\ ./classes/Inputs/TriaInput.cpp\ Index: /issm/trunk-jpl/src/c/analyses/bedslope_core.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/bedslope_core.cpp (revision 16381) +++ /issm/trunk-jpl/src/c/analyses/bedslope_core.cpp (revision 16382) @@ -14,7 +14,9 @@ /*parameters: */ bool save_results; + int meshtype; /*Recover some parameters: */ femmodel->parameters->FindParam(&save_results,SaveResultsEnum); + femmodel->parameters->FindParam(&meshtype,MeshTypeEnum); if(VerboseSolution()) _printf0_(" computing slope\n"); @@ -23,11 +25,13 @@ femmodel->SetCurrentConfiguration(BedSlopeAnalysisEnum,BedSlopeXAnalysisEnum); solutionsequence_linear(femmodel); - femmodel->SetCurrentConfiguration(BedSlopeAnalysisEnum,BedSlopeYAnalysisEnum); - solutionsequence_linear(femmodel); + if(meshtype!=Mesh2DverticalEnum){ + femmodel->SetCurrentConfiguration(BedSlopeAnalysisEnum,BedSlopeYAnalysisEnum); + solutionsequence_linear(femmodel); + } if(save_results){ if(VerboseSolution()) _printf0_(" saving results\n"); InputToResultx(femmodel,BedSlopeXEnum); - InputToResultx(femmodel,BedSlopeYEnum); + if(meshtype!=Mesh2DverticalEnum) InputToResultx(femmodel,BedSlopeYEnum); } Index: /issm/trunk-jpl/src/c/analyses/stressbalance_core.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/stressbalance_core.cpp (revision 16381) +++ /issm/trunk-jpl/src/c/analyses/stressbalance_core.cpp (revision 16382) @@ -49,5 +49,5 @@ /*Compute slopes: */ if(isSIA) surfaceslope_core(femmodel); - if(isFS && meshtype==Mesh3DEnum){ + if(isFS){ bedslope_core(femmodel); femmodel->SetCurrentConfiguration(StressbalanceAnalysisEnum); Index: /issm/trunk-jpl/src/c/classes/Elements/ElementHook.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/ElementHook.cpp (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Elements/ElementHook.cpp (revision 16382) @@ -128,2 +128,29 @@ } /*}}}*/ +/*FUNCTION ElementHook::SpawnSegHook{{{*/ +void ElementHook::SpawnSegHook(ElementHook* triahook,int index1,int index2){ + + triahook->numanalyses=this->numanalyses; + + int indices[2]; + indices[0]=index1; + indices[1]=index2; + + /*Spawn nodes hook*/ + triahook->hnodes=new Hook*[this->numanalyses]; + for(int i=0;inumanalyses;i++){ + /*Do not do anything if Hook is empty*/ + if (!this->hnodes[i] || this->hnodes[i]->GetNum()==0){ + triahook->hnodes[i]=NULL; + } + else{ + triahook->hnodes[i]=this->hnodes[i]->Spawn(indices,2); + } + } + + /*do not spawn hmaterial. material will be taken care of by Tria*/ + triahook->hmaterial=NULL; + triahook->hvertices=(Hook*)this->hvertices->Spawn(indices,2); + triahook->hmatpar=(Hook*)this->hmatpar->copy(); +} +/*}}}*/ Index: /issm/trunk-jpl/src/c/classes/Elements/ElementHook.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/ElementHook.h (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Elements/ElementHook.h (revision 16382) @@ -26,4 +26,5 @@ void SetHookNodes(int* node_ids,int numnodes,int analysis_counter); void SpawnTriaHook(ElementHook* triahook,int location); //3d only + void SpawnSegHook(ElementHook* triahook,int ndex1,int index2); //2d only void InitHookNeighbors(int* element_ids); //3d only }; Index: /issm/trunk-jpl/src/c/classes/Elements/Seg.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Seg.cpp (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Elements/Seg.cpp (revision 16382) @@ -139,2 +139,92 @@ } /*}}}*/ + +/*FUNCTION Seg::CreateMassMatrix {{{*/ +ElementMatrix* Seg::CreateMassMatrix(void){ + + /* Intermediaries */ + IssmDouble D,Jdet; + IssmDouble xyz_list[NUMVERTICES][2]; + + /*Fetch number of nodes and dof for this finite element*/ + int numnodes = this->NumberofNodes(); + + /*Initialize Element matrix and vectors*/ + ElementMatrix* Ke = new ElementMatrix(nodes,numnodes,this->parameters,NoneApproximationEnum); + IssmDouble* basis = xNew(numnodes); + + /*Retrieve all inputs and parameters*/ + GetVerticesCoordinates(&xyz_list[0][0],vertices,NUMVERTICES); + + /* Start looping on the number of gaussian points: */ + GaussSeg* gauss=new GaussSeg(2); + for(int ig=gauss->begin();igend();ig++){ + + gauss->GaussPoint(ig); + + GetNodalFunctions(basis,gauss); + GetJacobianDeterminant(&Jdet, &xyz_list[0][0],gauss); + D=gauss->weight*Jdet; + + TripleMultiply(basis,1,numnodes,1, + &D,1,1,0, + basis,1,numnodes,0, + &Ke->values[0],1); + } + + /*Clean up and return*/ + delete gauss; + xDelete(basis); + return Ke; +} +/*}}}*/ +/*FUNCTION Seg::CreatePVectorSlope {{{*/ +ElementVector* Seg::CreatePVectorSlope(void){ + + /*Intermediaries */ + int i,analysis_type; + IssmDouble Jdet; + IssmDouble xyz_list[NUMVERTICES][2]; + IssmDouble slope; + + /*Fetch number of nodes and dof for this finite element*/ + int numnodes = this->NumberofNodes(); + + /*Initialize Element vector*/ + ElementVector* pe = new ElementVector(nodes,numnodes,this->parameters); + IssmDouble* basis = xNew(numnodes); + + /*Retrieve all inputs and parameters*/ + parameters->FindParam(&analysis_type,AnalysisTypeEnum); + GetVerticesCoordinates(&xyz_list[0][0],vertices,NUMVERTICES); + Input* slope_input=NULL; + if(analysis_type==SurfaceSlopeXAnalysisEnum){ + slope_input=inputs->GetInput(SurfaceEnum); _assert_(slope_input); + } + else if(analysis_type==BedSlopeXAnalysisEnum){ + slope_input=inputs->GetInput(BedEnum); _assert_(slope_input); + } + else{ + _error_("Analysis "<begin();igend();ig++){ + + gauss->GaussPoint(ig); + + GetJacobianDeterminant(&Jdet, &xyz_list[0][0],gauss); + GetNodalFunctions(basis,gauss); + + slope_input->GetInputDerivativeValue(&slope,&xyz_list[0][0],gauss); + + for(i=0;ivalues[i]+=Jdet*gauss->weight*slope*basis[i]; + } + + /*Clean up and return*/ + xDelete(basis); + delete gauss; + return pe; +} +/*}}}*/ Index: /issm/trunk-jpl/src/c/classes/Elements/Seg.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Seg.h (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Elements/Seg.h (revision 16382) @@ -76,4 +76,5 @@ void CreateDVector(Vector* df){_error_("not implemented yet");}; void CreatePVector(Vector* pf){_error_("not implemented yet");}; + ElementVector* CreatePVectorSlope(void); void CreateJacobianMatrix(Matrix* Jff){_error_("not implemented yet");}; void Delta18oParameterization(void){_error_("not implemented yet");}; @@ -87,4 +88,19 @@ bool IsNodeOnShelfFromFlags(IssmDouble* flags){_error_("not implemented yet");}; bool NoIceInElement(){_error_("not implemented yet");}; + void GetInputListOnVertices(IssmDouble* pvalue,int enumtype){_error_("not implemented yet");}; + void GetInputListOnVertices(IssmDouble* pvalue,int enumtype,IssmDouble defaultvalue){_error_("not implemented yet");}; + void GetInputValue(IssmDouble* pvalue,Node* node,int enumtype){_error_("not implemented yet");}; + void GetMaterialInputValue(IssmDouble* pvalue,Node* node,int enumtype){_error_("not implemented yet");}; + #ifdef _HAVE_THERMAL_ + void UpdateBasalConstraintsEnthalpy(void){_error_("not implemented yet");}; + void ComputeBasalMeltingrate(void){_error_("not implemented yet");}; + void DrainWaterfraction(void){_error_("not implemented yet");}; + #endif + #ifdef _HAVE_HYDROLOGY_ + void GetHydrologyDCInefficientHmax(IssmDouble* ph_max, Node* innode){_error_("not implemented yet");}; + void GetHydrologyTransfer(Vector* transfer){_error_("not implemented yet");}; + void HydrologyEPLGetActive(Vector* active_vec){_error_("not implemented yet");}; + void HydrologyEPLGetMask(Vector* vec_mask){_error_("not implemented yet");}; + #endif void GetSolutionFromInputs(Vector* solution){_error_("not implemented yet");}; void GetVectorFromInputs(Vector* vector, int name_enum){_error_("not implemented yet");}; @@ -173,4 +189,6 @@ #endif /*}}}*/ + /*Seg specific routines:*/ + ElementMatrix* CreateMassMatrix(void); }; #endif /* _SEG_H */ Index: /issm/trunk-jpl/src/c/classes/Elements/SegRef.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/SegRef.cpp (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Elements/SegRef.cpp (revision 16382) @@ -48,5 +48,5 @@ /*Reference Element numerics*/ -/*FUNCTION SegRef::GetNodalFunctions{{{*/ +/*FUNCTION SegRef::GetNodalFunctions(IssmDouble* basis,GaussSeg* gauss){{{*/ void SegRef::GetNodalFunctions(IssmDouble* basis,GaussSeg* gauss){ /*This routine returns the values of the nodal functions at the gaussian point.*/ @@ -54,5 +54,14 @@ _assert_(basis); - switch(this->element_type){ + GetNodalFunctions(basis,gauss,this->element_type); +} +/*}}}*/ +/*FUNCTION SegRef::GetNodalFunctions(IssmDouble* basis,GaussSeg* gauss,int finiteelement){{{*/ +void SegRef::GetNodalFunctions(IssmDouble* basis,GaussSeg* gauss,int finiteelement){ + /*This routine returns the values of the nodal functions at the gaussian point.*/ + + _assert_(basis); + + switch(element_type){ case P1Enum: case P1DGEnum: basis[0]=(1.-gauss->coord1)/2.; @@ -60,6 +69,127 @@ return; default: + _error_("Element type "<element_type); + +} +/*}}}*/ +/*FUNCTION SegRef::GetNodalFunctionsDerivatives(IssmDouble* dbasis,IssmDouble* xyz_list, GaussSeg* gauss,int finiteelement){{{*/ +void SegRef::GetNodalFunctionsDerivatives(IssmDouble* dbasis,IssmDouble* xyz_list, GaussSeg* gauss,int finiteelement){ + + /*This routine returns the values of the nodal functions derivatives (with respect to the + * actual coordinate system): */ + IssmDouble Jinv; + + /*Fetch number of nodes for this finite element*/ + int numnodes = this->NumberofNodes(finiteelement); + + /*Get nodal functions derivatives in reference triangle*/ + IssmDouble* dbasis_ref=xNew(numnodes); + GetNodalFunctionsDerivativesReference(dbasis_ref,gauss,finiteelement); + + /*Get Jacobian invert: */ + GetJacobianInvert(&Jinv, xyz_list, gauss); + + /*Build dbasis: + * [dhi/dx]= Jinv*[dhi/dr] + */ + for(int i=0;i(dbasis_ref); + +} +/*}}}*/ +/*FUNCTION SegRef::GetNodalFunctionsDerivativesReference(IssmDouble* dbasis,GaussSeg* gauss){{{*/ +void SegRef::GetNodalFunctionsDerivativesReference(IssmDouble* dbasis,GaussSeg* gauss){ + /*This routine returns the values of the nodal functions derivatives (with respect to the + * natural coordinate system) at the gaussian point. */ + + GetNodalFunctionsDerivativesReference(dbasis,gauss,this->element_type); + +} +/*}}}*/ +/*FUNCTION SegRef::GetNodalFunctionsDerivativesReference(IssmDouble* dbasis,GaussSeg* gauss,int finiteelement){{{*/ +void SegRef::GetNodalFunctionsDerivativesReference(IssmDouble* dbasis,GaussSeg* gauss,int finiteelement){ + /*This routine returns the values of the nodal functions derivatives (with respect to the + * natural coordinate system) at the gaussian point. */ + + _assert_(dbasis && gauss); + + switch(finiteelement){ + case P1Enum: case P1DGEnum: + /*Nodal function 1*/ + dbasis[0] = -0.5; + /*Nodal function 2*/ + dbasis[1] = 0.5; + return; + default: _error_("Element type "<element_type)<<" not supported yet"); } + +} +/*}}}*/ +/*FUNCTION SegRef::GetInputDerivativeValue{{{*/ +void SegRef::GetInputDerivativeValue(IssmDouble* p, IssmDouble* plist,IssmDouble* xyz_list, GaussSeg* gauss){ + + /*From node values of parameter p (plist[0],plist[1]), return parameter derivative value at gaussian + * point specified by gauss_basis: + * dp/dx=plist[0]*dh1/dx+plist[1]*dh2/dx + * + * p is a vector already allocated. + */ + + /*Output*/ + IssmDouble dpx; + + /*Fetch number of nodes for this finite element*/ + int numnodes = this->NumberofNodes(); + + /*Get nodal functions derivatives*/ + IssmDouble* dbasis=xNew(1*numnodes); + GetNodalFunctionsDerivatives(dbasis,xyz_list,gauss); + + /*Calculate parameter for this Gauss point*/ + for(int i=0;i(dbasis); + *p=dpx; + +} +/*}}}*/ +/*FUNCTION SegRef::GetInputValue(IssmDouble* p, IssmDouble* plist, GaussSeg* gauss){{{*/ +void SegRef::GetInputValue(IssmDouble* p, IssmDouble* plist, GaussSeg* gauss){ + + GetInputValue(p,plist,gauss,this->element_type); +} +/*}}}*/ +/*FUNCTION SegRef::GetInputValue(IssmDouble* p, IssmDouble* plist, GaussSeg* gauss,int finiteelement){{{*/ +void SegRef::GetInputValue(IssmDouble* p, IssmDouble* plist, GaussSeg* gauss,int finiteelement){ + + /*Output*/ + IssmDouble value =0.; + + /*Fetch number of nodes for this finite element*/ + int numnodes = this->NumberofNodes(finiteelement); + + /*Get nodal functions*/ + IssmDouble* basis=xNew(numnodes); + GetNodalFunctions(basis, gauss,finiteelement); + + /*Calculate parameter for this Gauss point*/ + for(int i=0;i(basis); + *p = value; } /*}}}*/ @@ -103,2 +233,20 @@ } /*}}}*/ +/*FUNCTION SegRef::NumberofNodes(){{{*/ +int SegRef::NumberofNodes(void){ + + return this->NumberofNodes(this->element_type); +} +/*}}}*/ +/*FUNCTION SegRef::NumberofNodes(int finiteelement){{{*/ +int SegRef::NumberofNodes(int finiteelement){ + + switch(finiteelement){ + case P1Enum: return NUMNODESP1; + case P1DGEnum: return NUMNODESP1; + default: _error_("Element type "<element_type)<<" not supported yet"); + } + + return -1; +} +/*}}}*/ Index: /issm/trunk-jpl/src/c/classes/Elements/SegRef.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/SegRef.h (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Elements/SegRef.h (revision 16382) @@ -26,4 +26,15 @@ void GetJacobianInvert(IssmDouble* Jinv, IssmDouble* xyz_list,GaussSeg* gauss); void GetNodalFunctions(IssmDouble* basis,GaussSeg* gauss); + void GetNodalFunctions(IssmDouble* basis,GaussSeg* gauss,int finiteelement); + void GetNodalFunctionsDerivatives(IssmDouble* dbasis,IssmDouble* xyz_list, GaussSeg* gauss); + void GetNodalFunctionsDerivatives(IssmDouble* dbasis,IssmDouble* xyz_list, GaussSeg* gauss,int finiteelement); + void GetNodalFunctionsDerivativesReference(IssmDouble* dbasis,GaussSeg* gauss); + void GetNodalFunctionsDerivativesReference(IssmDouble* dbasis,GaussSeg* gauss,int finiteelement); + void GetInputDerivativeValue(IssmDouble* p, IssmDouble* plist,IssmDouble* xyz_list, GaussSeg* gauss); + void GetInputValue(IssmDouble* p, IssmDouble* plist, GaussSeg* gauss); + void GetInputValue(IssmDouble* p, IssmDouble* plist, GaussSeg* gauss,int finiteelement); + + int NumberofNodes(void); + int NumberofNodes(int finiteelement); }; #endif Index: /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp (revision 16382) @@ -246,5 +246,12 @@ #endif case BedSlopeXAnalysisEnum: case SurfaceSlopeXAnalysisEnum: case BedSlopeYAnalysisEnum: case SurfaceSlopeYAnalysisEnum: - return CreateMassMatrix(); + int meshtype; + parameters->FindParam(&meshtype,MeshTypeEnum); + if(meshtype==Mesh2DverticalEnum){ + return CreateBasalMassMatrix(); + } + else{ + return CreateMassMatrix(); + } break; #ifdef _HAVE_MASSTRANSPORT_ @@ -330,4 +337,20 @@ } /*}}}*/ +/*FUNCTION Tria::CreateBasalMassMatrix{{{*/ +ElementMatrix* Tria::CreateBasalMassMatrix(void){ + + if(!HasEdgeOnBed()) return NULL; + + int index1,index2; + this->EdgeOnBedIndices(&index1,&index2); + + Seg* seg=(Seg*)SpawnSeg(index1,index2); + ElementMatrix* Ke=seg->CreateMassMatrix(); + delete seg->material; delete seg; + + /*clean up and return*/ + return Ke; +} +/*}}}*/ /*FUNCTION Tria::CreateDVector {{{*/ void Tria::CreateDVector(Vector* df){ @@ -413,5 +436,12 @@ #endif case BedSlopeXAnalysisEnum: case SurfaceSlopeXAnalysisEnum: case BedSlopeYAnalysisEnum: case SurfaceSlopeYAnalysisEnum: - return CreatePVectorSlope(); + int meshtype; + parameters->FindParam(&meshtype,MeshTypeEnum); + if(meshtype==Mesh2DverticalEnum){ + return CreateBasalPVectorSlope(); + } + else{ + return CreatePVectorSlope(); + } break; #ifdef _HAVE_MASSTRANSPORT_ @@ -459,4 +489,20 @@ /*make compiler happy*/ return NULL; +} +/*}}}*/ +/*FUNCTION Tria::CreateBasalPVectorSlope{{{*/ +ElementVector* Tria::CreateBasalPVectorSlope(void){ + + if(!HasEdgeOnBed()) return NULL; + + int index1,index2; + this->EdgeOnBedIndices(&index1,&index2); + + Seg* seg=(Seg*)SpawnSeg(index1,index2); + ElementVector* pe=seg->CreatePVectorSlope(); + delete seg->material; delete seg; + + /*clean up and return*/ + return pe; } /*}}}*/ @@ -2459,4 +2505,34 @@ else this->nodes=NULL; +} +/*}}}*/ +/*FUNCTION Tria::SpawnSeg {{{*/ +Seg* Tria::SpawnSeg(int index1,int index2){ + + int analysis_counter; + + /*go into parameters and get the analysis_counter: */ + this->parameters->FindParam(&analysis_counter,AnalysisCounterEnum); + + /*Create Seg*/ + Seg* seg=new Seg(); + seg->id=this->id; + seg->inputs=(Inputs*)this->inputs->SpawnSegInputs(index1,index2); + seg->parameters=this->parameters; + seg->element_type=P1Enum; //Only P1 CG for now (TO BE CHANGED) + this->SpawnSegHook(dynamic_cast(seg),index1,index2); + + /*Spawn material*/ + seg->material=(Material*)this->material->copy(); + delete seg->material->inputs; + seg->material->inputs=(Inputs*)this->material->inputs->SpawnSegInputs(index1,index2); + + /*recover nodes, material and matpar: */ + seg->nodes = (Node**)seg->hnodes[analysis_counter]->deliverp(); + seg->vertices = (Vertex**)seg->hvertices->deliverp(); + seg->matpar = (Matpar*)seg->hmatpar->delivers(); + + /*Return new Seg*/ + return seg; } /*}}}*/ Index: /issm/trunk-jpl/src/c/classes/Elements/Tria.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Tria.h (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Elements/Tria.h (revision 16382) @@ -18,4 +18,5 @@ class Material; class Matpar; +class Seg; class ElementMatrix; class ElementVector; @@ -190,4 +191,5 @@ ElementMatrix* CreateKMatrixFreeSurfaceBase(void); ElementMatrix* CreateMassMatrix(void); + ElementMatrix* CreateBasalMassMatrix(void); ElementVector* CreatePVector(void); ElementVector* CreatePVectorBalancethickness(void); @@ -203,4 +205,5 @@ ElementVector* CreatePVectorFreeSurfaceBase(void); ElementVector* CreatePVectorSlope(void); + ElementVector* CreateBasalPVectorSlope(void); IssmDouble GetArea(void); void GetAreaCoordinates(IssmDouble *area_coordinates,IssmDouble xyz_zero[3][3],IssmDouble xyz_list[3][3],int numpoints); @@ -228,4 +231,5 @@ bool IsInput(int name); void SetClone(int* minranks); + Seg* SpawnSeg(int index1,int index2); void SurfaceNormal(IssmDouble* surface_normal, IssmDouble xyz_list[3][3]); Index: /issm/trunk-jpl/src/c/classes/Inputs/BoolInput.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/BoolInput.cpp (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Inputs/BoolInput.cpp (revision 16382) @@ -73,4 +73,19 @@ /*FUNCTION BoolInput::SpawnTriaInput{{{*/ Input* BoolInput::SpawnTriaInput(int location){ + + /*output*/ + BoolInput* outinput=new BoolInput(); + + /*only copy current value*/ + outinput->enum_type=this->enum_type; + outinput->value=this->value; + + /*Assign output*/ + return outinput; + +} +/*}}}*/ +/*FUNCTION BoolInput::SpawnSegInput{{{*/ +Input* BoolInput::SpawnSegInput(int index1,int index2){ /*output*/ Index: /issm/trunk-jpl/src/c/classes/Inputs/BoolInput.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/BoolInput.h (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Inputs/BoolInput.h (revision 16382) @@ -35,4 +35,5 @@ int InstanceEnum(); Input* SpawnTriaInput(int location); + Input* SpawnSegInput(int index1,int index2); Input* PointwiseDivide(Input* inputB){_error_("not implemented yet");}; Input* PointwiseMin(Input* inputB){_error_("not implemented yet");}; @@ -52,4 +53,5 @@ void GetInputValue(IssmDouble* pvalue,GaussTria* gauss ,int index){_error_("not implemented yet");}; void GetInputValue(IssmDouble* pvalue,GaussPenta* gauss ,int index){_error_("not implemented yet");}; + void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussSeg* gauss){_error_("not implemented yet");}; void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussTria* gauss); void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussPenta* gauss); Index: /issm/trunk-jpl/src/c/classes/Inputs/ControlInput.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/ControlInput.cpp (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Inputs/ControlInput.cpp (revision 16382) @@ -186,4 +186,8 @@ return values->SpawnTriaInput(location); }/*}}}*/ +/*FUNCTION ControlInput::SpawnSegInput{{{*/ +Input* ControlInput::SpawnSegInput(int index1,int index2){ + return values->SpawnSegInput(index1,index2); +}/*}}}*/ /*FUNCTION ControlInput::SpawnGradient{{{*/ ElementResult* ControlInput::SpawnGradient(int step, IssmDouble time){ Index: /issm/trunk-jpl/src/c/classes/Inputs/ControlInput.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/ControlInput.h (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Inputs/ControlInput.h (revision 16382) @@ -39,4 +39,5 @@ int InstanceEnum(); Input* SpawnTriaInput(int location); + Input* SpawnSegInput(int index1,int index2); Input* PointwiseDivide(Input* inputB){_error_("not implemented yet");}; Input* PointwiseMin(Input* inputB){_error_("not implemented yet");}; @@ -57,4 +58,5 @@ void GetInputValue(IssmDouble* pvalue,GaussTria* gauss ,int index){_error_("not implemented yet");}; void GetInputValue(IssmDouble* pvalue,GaussPenta* gauss ,int index){_error_("not implemented yet");}; + void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussSeg* gauss){_error_("not implemented yet");}; void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussTria* gauss); void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussPenta* gauss); Index: /issm/trunk-jpl/src/c/classes/Inputs/DatasetInput.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/DatasetInput.cpp (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Inputs/DatasetInput.cpp (revision 16382) @@ -96,4 +96,22 @@ } /*}}}*/ +/*FUNCTION DatasetInput::SpawnSegInput{{{*/ +Input* DatasetInput::SpawnSegInput(int index1,int index2){ + + /*output*/ + DatasetInput* outinput=NULL; + + /*Create new Datasetinput (copy of current input)*/ + outinput=new DatasetInput(); + outinput->enum_type=this->enum_type; + outinput->inputs=dynamic_cast(this->inputs->SpawnSegInputs(index1,index2)); + outinput->numids=this->numids; + outinput->ids=xNew(this->numids); + xMemCpy(outinput->ids,this->ids,this->numids); + + /*Assign output*/ + return outinput; +} +/*}}}*/ /*DatasetInput management*/ Index: /issm/trunk-jpl/src/c/classes/Inputs/DatasetInput.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/DatasetInput.h (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Inputs/DatasetInput.h (revision 16382) @@ -37,4 +37,5 @@ int InstanceEnum(); Input* SpawnTriaInput(int location); + Input* SpawnSegInput(int index1,int index2); Input* PointwiseDivide(Input* inputB){_error_("not implemented yet");}; Input* PointwiseMin(Input* inputB){_error_("not implemented yet");}; @@ -54,4 +55,5 @@ void GetInputValue(IssmDouble* pvalue,GaussTria* gauss ,int index); void GetInputValue(IssmDouble* pvalue,GaussPenta* gauss ,int index); + void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussSeg* gauss){_error_("not implemented yet");}; void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussTria* gauss){_error_("not implemented yet");}; void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussPenta* gauss){_error_("not implemented yet");}; Index: /issm/trunk-jpl/src/c/classes/Inputs/DoubleInput.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/DoubleInput.cpp (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Inputs/DoubleInput.cpp (revision 16382) @@ -86,4 +86,19 @@ } /*}}}*/ +/*FUNCTION DoubleInput::SpawnSegInput{{{*/ +Input* DoubleInput::SpawnSegInput(int index1,int index2){ + + /*output*/ + DoubleInput* outinput=new DoubleInput(); + + /*only copy current value*/ + outinput->enum_type=this->enum_type; + outinput->value=this->value; + + /*Assign output*/ + return outinput; + +} +/*}}}*/ /*FUNCTION DoubleInput::SpawnResult{{{*/ ElementResult* DoubleInput::SpawnResult(int step, IssmDouble time){ Index: /issm/trunk-jpl/src/c/classes/Inputs/DoubleInput.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/DoubleInput.h (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Inputs/DoubleInput.h (revision 16382) @@ -34,4 +34,5 @@ int InstanceEnum(); Input* SpawnTriaInput(int location); + Input* SpawnSegInput(int index1,int index2); Input* PointwiseDivide(Input* inputB); Input* PointwiseMin(Input* inputB); @@ -51,4 +52,5 @@ void GetInputValue(IssmDouble* pvalue,GaussTria* gauss ,int index){_error_("not implemented yet");}; void GetInputValue(IssmDouble* pvalue,GaussPenta* gauss ,int index){_error_("not implemented yet");}; + void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussSeg* gauss){_error_("not implemented yet");}; void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussTria* gauss); void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussPenta* gauss); Index: /issm/trunk-jpl/src/c/classes/Inputs/Input.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/Input.h (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Inputs/Input.h (revision 16382) @@ -13,4 +13,5 @@ class ElementResult; class GaussTria; +class GaussSeg; class Parameters; class GaussPenta; @@ -34,4 +35,5 @@ virtual void GetInputValue(IssmDouble* pvalue,GaussTria* gauss ,int index)=0; virtual void GetInputValue(IssmDouble* pvalue,GaussPenta* gauss,int index)=0; + virtual void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussSeg* gauss)=0; virtual void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussTria* gauss)=0; virtual void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussPenta* gauss)=0; @@ -65,4 +67,5 @@ virtual Input* SpawnTriaInput(int location)=0; + virtual Input* SpawnSegInput(int index1,int index2)=0; virtual Input* PointwiseDivide(Input* inputB)=0; virtual Input* PointwiseMax(Input* inputmax)=0; Index: /issm/trunk-jpl/src/c/classes/Inputs/Inputs.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/Inputs.cpp (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Inputs/Inputs.cpp (revision 16382) @@ -395,4 +395,30 @@ } /*}}}*/ +/*FUNCTION Inputs::SpawnSegInputs{{{*/ +Inputs* Inputs::SpawnSegInputs(int index1,int index2){ + + /*Intermediary*/ + vector::iterator object; + Input* inputin=NULL; + Input* inputout=NULL; + + /*Output*/ + Inputs* newinputs=new Inputs(); + + /*Go through inputs and call Spawn function*/ + for ( object=objects.begin() ; object < objects.end(); object++ ){ + + /*Create new input*/ + inputin=dynamic_cast(*object); + inputout=inputin->SpawnSegInput(index1,index2); + + /*Add input to new inputs*/ + newinputs->AddObject(inputout); + } + + /*Assign output pointer*/ + return newinputs; +} +/*}}}*/ /*FUNCTION Inputs::AXPY{{{*/ void Inputs::AXPY(int inputy_enum, IssmDouble scalar, int inputx_enum){ Index: /issm/trunk-jpl/src/c/classes/Inputs/Inputs.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/Inputs.h (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Inputs/Inputs.h (revision 16382) @@ -29,4 +29,6 @@ Input* GetInput(int enum_name); Inputs* SpawnTriaInputs(int position); + Inputs* SpawnSegInputs(int index1,int index2); + Inputs* SpawnSegInputs(int position); void AXPY(int inputy_enum, IssmDouble scalar, int inputx_enum); IssmDouble InfinityNorm(int enumtype); Index: /issm/trunk-jpl/src/c/classes/Inputs/IntInput.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/IntInput.cpp (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Inputs/IntInput.cpp (revision 16382) @@ -73,4 +73,18 @@ /*FUNCTION IntInput::SpawnTriaInput{{{*/ Input* IntInput::SpawnTriaInput(int location){ + + /*output*/ + IntInput* outinput=new IntInput(); + + /*only copy current value*/ + outinput->enum_type=this->enum_type; + outinput->value=this->value; + + /*Assign output*/ + return outinput; +} +/*}}}*/ +/*FUNCTION IntInput::SpawnSegInput{{{*/ +Input* IntInput::SpawnSegInput(int index1,int index2){ /*output*/ Index: /issm/trunk-jpl/src/c/classes/Inputs/IntInput.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/IntInput.h (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Inputs/IntInput.h (revision 16382) @@ -35,4 +35,5 @@ int InstanceEnum(); Input* SpawnTriaInput(int location); + Input* SpawnSegInput(int index1,int index2); Input* PointwiseDivide(Input* inputB){_error_("not implemented yet");}; Input* PointwiseMin(Input* inputB){_error_("not implemented yet");}; @@ -52,4 +53,5 @@ void GetInputValue(IssmDouble* pvalue,GaussTria* gauss ,int index){_error_("not implemented yet");}; void GetInputValue(IssmDouble* pvalue,GaussPenta* gauss ,int index){_error_("not implemented yet");}; + void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussSeg* gauss){_error_("not implemented yet");}; void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussTria* gauss); void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussPenta* gauss); Index: /issm/trunk-jpl/src/c/classes/Inputs/PentaInput.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/PentaInput.cpp (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Inputs/PentaInput.cpp (revision 16382) @@ -124,4 +124,10 @@ } /*}}}*/ +/*FUNCTION PentaInput::SpawnSegInput{{{*/ +Input* PentaInput::SpawnSegInput(int index1,int index2){ + + _error_("not supported"); +} +/*}}}*/ /*FUNCTION PentaInput::SpawnResult{{{*/ ElementResult* PentaInput::SpawnResult(int step, IssmDouble time){ Index: /issm/trunk-jpl/src/c/classes/Inputs/PentaInput.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/PentaInput.h (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Inputs/PentaInput.h (revision 16382) @@ -35,4 +35,5 @@ int InstanceEnum(); Input* SpawnTriaInput(int location); + Input* SpawnSegInput(int index1,int index2); Input* PointwiseDivide(Input* inputB); Input* PointwiseMin(Input* inputB); @@ -51,4 +52,5 @@ void GetInputValue(IssmDouble* pvalue,GaussTria* gauss ,int index){_error_("not implemented yet");}; void GetInputValue(IssmDouble* pvalue,GaussPenta* gauss ,int index){_error_("not implemented yet");}; + void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussSeg* gauss){_error_("not implemented yet");}; void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussTria* gauss){_error_("not implemented yet");}; void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussPenta* gauss); Index: /issm/trunk-jpl/src/c/classes/Inputs/SegInput.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/SegInput.cpp (revision 16382) +++ /issm/trunk-jpl/src/c/classes/Inputs/SegInput.cpp (revision 16382) @@ -0,0 +1,106 @@ +/*!\file SegInput.c + * \brief: implementation of the SegInput object + */ + +#ifdef HAVE_CONFIG_H + #include +#else +#error "Cannot compile with HAVE_CONFIG_H symbol! run configure first!" +#endif + +#include "../classes.h" +#include "../../shared/shared.h" + +/*SegInput constructors and destructor*/ +/*FUNCTION SegInput::SegInput(){{{*/ +SegInput::SegInput(){ + values = NULL; +} +/*}}}*/ +/*FUNCTION SegInput::SegInput(int in_enum_type,IssmDouble* invalues,element_type_in){{{*/ +SegInput::SegInput(int in_enum_type,IssmDouble* in_values,int element_type_in) + :SegRef(1) +{ + + /*Set SegRef*/ + this->SetElementType(element_type_in,0); + this->element_type=element_type_in; + + /*Set Enum*/ + enum_type=in_enum_type; + + /*Set values*/ + this->values=xNew(this->NumberofNodes()); + for(int i=0;iNumberofNodes();i++) values[i]=in_values[i]; +} +/*}}}*/ +/*FUNCTION SegInput::~SegInput(){{{*/ +SegInput::~SegInput(){ + xDelete(this->values); +} +/*}}}*/ + +/*Object virtual functions definitions:*/ +/*FUNCTION SegInput::Echo {{{*/ +void SegInput::Echo(void){ + this->DeepEcho(); +} +/*}}}*/ +/*FUNCTION SegInput::DeepEcho{{{*/ +void SegInput::DeepEcho(void){ + + _printf_("SegInput:\n"); + _printf_(" enum: " << this->enum_type << " (" << EnumToStringx(this->enum_type) << ")\n"); + _printf_(" values: ["); + for(int i=0;iNumberofNodes();i++) _printf_(" "<values[i]); + _printf_("]\n"); +} +/*}}}*/ +/*FUNCTION SegInput::Id{{{*/ +int SegInput::Id(void){ return -1; } +/*}}}*/ +/*FUNCTION SegInput::ObjectEnum{{{*/ +int SegInput::ObjectEnum(void){ + + return SegInputEnum; + +} +/*}}}*/ +/*FUNCTION SegInput::copy{{{*/ +Object* SegInput::copy() { + + return new SegInput(this->enum_type,this->values,this->element_type); + +} +/*}}}*/ + +/*SegInput management*/ +/*FUNCTION SegInput::InstanceEnum{{{*/ +int SegInput::InstanceEnum(void){ + + return this->enum_type; + +} +/*}}}*/ + +/*Object functions*/ +/*FUNCTION SegInput::GetInputValue(IssmDouble* pvalue,GaussSeg* gauss){{{*/ +void SegInput::GetInputValue(IssmDouble* pvalue,GaussSeg* gauss){ + + /*Call SegRef function*/ + SegRef::GetInputValue(pvalue,&values[0],gauss); + +} +/*}}}*/ +/*FUNCTION SegInput::GetInputDerivativeValue(IssmDouble* p, IssmDouble* xyz_list, GaussSeg* gauss){{{*/ +void SegInput::GetInputDerivativeValue(IssmDouble* p, IssmDouble* xyz_list, GaussSeg* gauss){ + + /*Call SegRef function*/ + SegRef::GetInputDerivativeValue(p,&values[0],xyz_list,gauss); +} +/*}}}*/ +/*FUNCTION SegInput::Configure{{{*/ +void SegInput::Configure(Parameters* parameters){ + /*do nothing: */ +} +/*}}}*/ Index: /issm/trunk-jpl/src/c/classes/Inputs/SegInput.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/SegInput.h (revision 16382) +++ /issm/trunk-jpl/src/c/classes/Inputs/SegInput.h (revision 16382) @@ -0,0 +1,87 @@ +/*! \file SegInput.h + * \brief: header file for SegInput object + */ + +#ifndef _SEGINPUT_H_ +#define _SEGINPUT_H_ + +/*Headers:*/ +/*{{{*/ +#include "./Input.h" +#include "../Elements/SegRef.h" +class GaussSeg; +class GaussPenta; +/*}}}*/ + +class SegInput: public Input,public SegRef{ + + public: + int enum_type; + IssmDouble* values; + + /*SegInput constructors, destructors*/ + SegInput(); + SegInput(int enum_type,IssmDouble* values,int element_type_in); + ~SegInput(); + + /*Object virtual functions definitions*/ + void Echo(); + void DeepEcho(); + int Id(); + int ObjectEnum(); + Object *copy(); + + /*SegInput management:*/ + int InstanceEnum(); + Input* SpawnTriaInput(int location){_error_("not supported yet");}; + Input* SpawnSegInput(int index1,int index2){_error_("not implemented yet");}; + Input* PointwiseDivide(Input* inputB){_error_("not supported yet");}; + Input* PointwiseMin(Input* inputB){_error_("not supported yet");}; + Input* PointwiseMax(Input* inputB){_error_("not supported yet");}; + ElementResult* SpawnResult(int step, IssmDouble time){_error_("not supported yet");}; + void AddTimeValues(IssmDouble* values,int step,IssmDouble time){_error_("not supported yet");}; + void Configure(Parameters* parameters); + + /*numerics*/ + void GetInputValue(bool* pvalue){_error_("not implemented yet");} + void GetInputValue(int* pvalue){_error_("not implemented yet");} + void GetInputValue(IssmDouble* pvalue){_error_("not implemented yet");} + void GetInputValue(IssmDouble* pvalue,GaussSeg* gauss); + void GetInputValue(IssmDouble* pvalue,GaussTria* gauss){_error_("not implemented yet");}; + void GetInputValue(IssmDouble* pvalue,GaussPenta* gauss){_error_("not implemented yet");}; + void GetInputValue(IssmDouble* pvalue,GaussTria* gauss,IssmDouble time){_error_("not implemented yet");}; + void GetInputValue(IssmDouble* pvalue,GaussPenta* gauss,IssmDouble time){_error_("not implemented yet");}; + void GetInputValue(IssmDouble* pvalue,GaussTria* gauss ,int index){_error_("not implemented yet");}; + void GetInputValue(IssmDouble* pvalue,GaussPenta* gauss,int index){_error_("not implemented yet");}; + void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussSeg* gauss); + void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussTria* gauss){_error_("not implemented yet");}; + void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussPenta* gauss){_error_("not implemented yet");}; + void GetInputAverage(IssmDouble* pvalue){_error_("not implemented yet");}; + void GetInputAllTimeAverages(IssmDouble** pvalues,IssmDouble** ptimes, int* pnumtimes){_error_("not implemented yet");}; + void GetInputUpToCurrentTimeAverages(IssmDouble** pvalues, IssmDouble** ptimes, int* pnumtimes, IssmDouble currenttime){_error_("not implemented yet");}; + void GetVxStrainRate2d(IssmDouble* epsilonvx,IssmDouble* xyz_list, GaussTria* gauss){_error_("not implemented yet");}; + void GetVyStrainRate2d(IssmDouble* epsilonvy,IssmDouble* xyz_list, GaussTria* gauss){_error_("not implemented yet");}; + void GetVxStrainRate3d(IssmDouble* epsilonvx,IssmDouble* xyz_list, GaussPenta* gauss){_error_("not implemented yet");}; + void GetVyStrainRate3d(IssmDouble* epsilonvy,IssmDouble* xyz_list, GaussPenta* gauss){_error_("not implemented yet");}; + void GetVzStrainRate3d(IssmDouble* epsilonvz,IssmDouble* xyz_list, GaussPenta* gauss){_error_("not implemented yet");}; + void GetVxStrainRate3dHO(IssmDouble* epsilonvx,IssmDouble* xyz_list, GaussPenta* gauss){_error_("not implemented yet");}; + void GetVyStrainRate3dHO(IssmDouble* epsilonvy,IssmDouble* xyz_list, GaussPenta* gauss){_error_("not implemented yet");}; + void ChangeEnum(int newenumtype){_error_("not implemented yet");}; + + void SquareMin(IssmDouble* psquaremin,Parameters* parameters){_error_("not implemented yet");}; + void ConstrainMin(IssmDouble minimum){_error_("not implemented yet");}; + void Set(IssmDouble setvalue){_error_("not implemented yet");}; + void Scale(IssmDouble scale_factor){_error_("not implemented yet");}; + void AXPY(Input* xinput,IssmDouble scalar){_error_("not implemented yet");}; + void Constrain(IssmDouble cm_min, IssmDouble cm_max){_error_("not implemented yet");}; + IssmDouble InfinityNorm(void){_error_("not implemented yet");}; + IssmDouble Max(void){_error_("not implemented yet");}; + IssmDouble MaxAbs(void){_error_("not implemented yet");}; + IssmDouble Min(void){_error_("not implemented yet");}; + IssmDouble MinAbs(void){_error_("not implemented yet");}; + void Extrude(void){_error_("not supported yet");}; + void VerticallyIntegrate(Input* thickness_input){_error_("not supported yet");}; + void GetVectorFromInputs(Vector* vector,int* doflist){_error_("not implemented yet");}; + +}; +#endif /* _SEGINPUT_H */ Index: /issm/trunk-jpl/src/c/classes/Inputs/TransientInput.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/TransientInput.cpp (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Inputs/TransientInput.cpp (revision 16382) @@ -120,4 +120,24 @@ xMemCpy(outinput->timesteps,this->timesteps,this->numtimesteps); outinput->inputs=(Inputs*)this->inputs->SpawnTriaInputs(location); + outinput->parameters=this->parameters; + + /*Assign output*/ + return outinput; + +} +/*}}}*/ +/*FUNCTION TransientInput::SpawnSegInput{{{*/ +Input* TransientInput::SpawnSegInput(int index1,int index2){ + + /*output*/ + TransientInput* outinput=NULL; + + /*Create new Transientinput (copy of current input)*/ + outinput=new TransientInput(); + outinput->enum_type=this->enum_type; + outinput->numtimesteps=this->numtimesteps; + outinput->timesteps=xNew(this->numtimesteps); + xMemCpy(outinput->timesteps,this->timesteps,this->numtimesteps); + outinput->inputs=(Inputs*)this->inputs->SpawnSegInputs(index1,index2); outinput->parameters=this->parameters; Index: /issm/trunk-jpl/src/c/classes/Inputs/TransientInput.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/TransientInput.h (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Inputs/TransientInput.h (revision 16382) @@ -39,4 +39,5 @@ int InstanceEnum(); Input* SpawnTriaInput(int location); + Input* SpawnSegInput(int index1,int index2); Input* PointwiseDivide(Input* forcingB){_error_("not implemented yet");}; Input* PointwiseMin(Input* forcingB){_error_("not implemented yet");}; @@ -55,4 +56,5 @@ void GetInputValue(IssmDouble* pvalue,GaussTria* gauss ,int index){_error_("not implemented yet");}; void GetInputValue(IssmDouble* pvalue,GaussPenta* gauss ,int index){_error_("not implemented yet");}; + void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussSeg* gauss){_error_("not implemented yet");}; void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussTria* gauss); void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussPenta* gauss){_error_("not implemented yet");}; Index: /issm/trunk-jpl/src/c/classes/Inputs/TriaInput.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/TriaInput.cpp (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Inputs/TriaInput.cpp (revision 16382) @@ -97,4 +97,24 @@ } /*}}}*/ +/*FUNCTION SegInput::SpawnSegInput{{{*/ +Input* TriaInput::SpawnSegInput(int index1,int index2){ + + /*output*/ + SegInput* outinput=NULL; + IssmDouble newvalues[2]; //Assume P1 interpolation only for now + + /*Create arrow of indices depending on location (0=base 1=surface)*/ + + newvalues[0]=this->values[index1]; + newvalues[1]=this->values[index2]; + + /*Create new Seg input*/ + outinput=new SegInput(this->enum_type,&newvalues[0],P1Enum); + + /*Assign output*/ + return outinput; + +} +/*}}}*/ /*FUNCTION TriaInput::SpawnResult{{{*/ ElementResult* TriaInput::SpawnResult(int step, IssmDouble time){ Index: /issm/trunk-jpl/src/c/classes/Inputs/TriaInput.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/TriaInput.h (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Inputs/TriaInput.h (revision 16382) @@ -35,4 +35,5 @@ int InstanceEnum(); Input* SpawnTriaInput(int location); + Input* SpawnSegInput(int index1,int index2); Input* PointwiseDivide(Input* inputB); Input* PointwiseMin(Input* inputB); @@ -52,4 +53,5 @@ void GetInputValue(IssmDouble* pvalue,GaussTria* gauss ,int index){_error_("not implemented yet");}; void GetInputValue(IssmDouble* pvalue,GaussPenta* gauss,int index){_error_("not implemented yet");}; + void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussSeg* gauss){_error_("not implemented yet");}; void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussTria* gauss); void GetInputDerivativeValue(IssmDouble* derivativevalues, IssmDouble* xyz_list, GaussPenta* gauss){_error_("not implemented yet");}; Index: /issm/trunk-jpl/src/c/classes/Node.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Node.cpp (revision 16381) +++ /issm/trunk-jpl/src/c/classes/Node.cpp (revision 16382) @@ -95,5 +95,5 @@ analysis_type==HydrologyDCEfficientAnalysisEnum ){ - if(iomodel->meshtype==Mesh3DEnum){ + if(iomodel->meshtype==Mesh3DEnum || iomodel->meshtype==Mesh2DverticalEnum){ /*On a 3d mesh, we may have collapsed elements, hence dead nodes. Freeze them out: */ _assert_(iomodel->Data(MeshVertexonbedEnum)); Index: /issm/trunk-jpl/src/c/classes/classes.h =================================================================== --- /issm/trunk-jpl/src/c/classes/classes.h (revision 16381) +++ /issm/trunk-jpl/src/c/classes/classes.h (revision 16382) @@ -59,4 +59,5 @@ #include "./Inputs/PentaInput.h" #include "./Inputs/TriaInput.h" +#include "./Inputs/SegInput.h" #include "./Inputs/ControlInput.h" #include "./Inputs/DatasetInput.h" Index: /issm/trunk-jpl/src/c/modules/ModelProcessorx/BedSlope/CreateNodesBedSlope.cpp =================================================================== --- /issm/trunk-jpl/src/c/modules/ModelProcessorx/BedSlope/CreateNodesBedSlope.cpp (revision 16381) +++ /issm/trunk-jpl/src/c/modules/ModelProcessorx/BedSlope/CreateNodesBedSlope.cpp (revision 16382) @@ -11,5 +11,10 @@ void CreateNodesBedSlope(Nodes** pnodes, IoModel* iomodel){ - if(iomodel->meshtype==Mesh3DEnum) iomodel->FetchData(2,MeshVertexonbedEnum,MeshVertexonsurfaceEnum); + if(iomodel->meshtype==Mesh3DEnum){ + iomodel->FetchData(2,MeshVertexonbedEnum,MeshVertexonsurfaceEnum); + } + else if(iomodel->meshtype==Mesh2DverticalEnum){ + iomodel->FetchData(1,MeshVertexonbedEnum); + } CreateNodes(pnodes,iomodel,BedSlopeAnalysisEnum,P1Enum); iomodel->DeleteData(2,MeshVertexonbedEnum,MeshVertexonsurfaceEnum); Index: /issm/trunk-jpl/src/c/modules/ModelProcessorx/BedSlope/UpdateElementsBedSlope.cpp =================================================================== --- /issm/trunk-jpl/src/c/modules/ModelProcessorx/BedSlope/UpdateElementsBedSlope.cpp (revision 16381) +++ /issm/trunk-jpl/src/c/modules/ModelProcessorx/BedSlope/UpdateElementsBedSlope.cpp (revision 16382) @@ -28,3 +28,6 @@ iomodel->FetchDataToInput(elements,MeshElementonsurfaceEnum); } + if(iomodel->meshtype==Mesh2DverticalEnum){ + iomodel->FetchDataToInput(elements,MeshVertexonbedEnum); + } } Index: /issm/trunk-jpl/src/c/modules/ModelProcessorx/SurfaceSlope/CreateNodesSurfaceSlope.cpp =================================================================== --- /issm/trunk-jpl/src/c/modules/ModelProcessorx/SurfaceSlope/CreateNodesSurfaceSlope.cpp (revision 16381) +++ /issm/trunk-jpl/src/c/modules/ModelProcessorx/SurfaceSlope/CreateNodesSurfaceSlope.cpp (revision 16382) @@ -11,5 +11,10 @@ void CreateNodesSurfaceSlope(Nodes** pnodes, IoModel* iomodel){ - if(iomodel->meshtype==Mesh3DEnum) iomodel->FetchData(2,MeshVertexonbedEnum,MeshVertexonsurfaceEnum); + if(iomodel->meshtype==Mesh3DEnum){ + iomodel->FetchData(2,MeshVertexonbedEnum,MeshVertexonsurfaceEnum); + } + else if(iomodel->meshtype==Mesh2DverticalEnum){ + iomodel->FetchData(1,MeshVertexonbedEnum); + } CreateNodes(pnodes,iomodel,SurfaceSlopeAnalysisEnum,P1Enum); iomodel->DeleteData(2,MeshVertexonbedEnum,MeshVertexonsurfaceEnum);