Index: /issm/trunk/src/c/EnumDefinitions/EnumDefinitions.h =================================================================== --- /issm/trunk/src/c/EnumDefinitions/EnumDefinitions.h (revision 8802) +++ /issm/trunk/src/c/EnumDefinitions/EnumDefinitions.h (revision 8803) @@ -375,5 +375,4 @@ IsMacAyealPattynEnum, IsStokesEnum, - KffEnum, LatentHeatEnum, LowmemEnum, Index: /issm/trunk/src/c/Makefile.am =================================================================== --- /issm/trunk/src/c/Makefile.am (revision 8802) +++ /issm/trunk/src/c/Makefile.am (revision 8803) @@ -273,6 +273,4 @@ ./objects/Params/StringParam.cpp\ ./objects/Params/StringParam.h\ - ./objects/NodeSets.cpp\ - ./objects/NodeSets.h\ ./Container/Container.h\ ./Container/Constraints.h\ @@ -406,8 +404,6 @@ ./io/WriteDataToDisk.cpp\ ./io/IoModelFetchData.cpp\ - ./io/WriteNodeSets.cpp\ ./io/WriteParams.cpp\ ./io/FetchParams.cpp\ - ./io/FetchNodeSets.cpp\ ./io/OptionParse.cpp\ ./io/MatlabNArrayToNArray.cpp\ @@ -564,4 +560,6 @@ ./modules/Orthx/Orthx.h\ ./modules/Orthx/Orthx.cpp\ + ./modules/UpdateSpcsx/UpdateSpcsx.h\ + ./modules/UpdateSpcsx/UpdateSpcsx.cpp\ ./modules/Gradjx/Gradjx.h\ ./modules/Gradjx/Gradjx.cpp\ @@ -594,7 +592,4 @@ ./modules/ComputeStrainRatex/ComputeStrainRatex.h\ ./modules/ComputeStrainRatex/ComputeStrainRatex.cpp\ - ./modules/BuildNodeSetsx/BuildNodeSetsx.h\ - ./modules/BuildNodeSetsx/BuildNodeSetsx.cpp\ - ./modules/BuildNodeSetsx/PartitionSets.cpp\ ./modules/SpcNodesx/SpcNodesx.h\ ./modules/SpcNodesx/SpcNodesx.cpp\ @@ -634,8 +629,4 @@ ./modules/Reducevectorgtofx/Reducevectorgtofx.cpp\ ./modules/Reducevectorgtofx/Reducevectorgtofx.h\ - ./modules/Reducematrixfromgtofx/Reducematrixfromgtofx.cpp\ - ./modules/Reducematrixfromgtofx/Reducematrixfromgtofx.h\ - ./modules/Reduceloadfromgtofx/Reduceloadfromgtofx.h\ - ./modules/Reduceloadfromgtofx/Reduceloadfromgtofx.cpp\ ./modules/Reduceloadx/Reduceloadx.h\ ./modules/Reduceloadx/Reduceloadx.cpp\ @@ -670,4 +661,6 @@ ./modules/Solverx/Solverx.cpp\ ./modules/Solverx/Solverx.h\ + ./modules/VecMergex/VecMergex.cpp\ + ./modules/VecMergex/VecMergex.h\ ./modules/Mergesolutionfromftogx/Mergesolutionfromftogx.cpp\ ./modules/Mergesolutionfromftogx/Mergesolutionfromftogx.h\ @@ -957,6 +950,4 @@ ./objects/Params/StringParam.cpp\ ./objects/Params/StringParam.h\ - ./objects/NodeSets.cpp\ - ./objects/NodeSets.h\ ./Container/Container.h\ ./Container/Constraints.h\ @@ -1088,7 +1079,5 @@ ./io/WriteDataToDisk.cpp\ ./io/IoModelFetchData.cpp\ - ./io/WriteNodeSets.cpp\ ./io/WriteParams.cpp\ - ./io/FetchNodeSets.cpp\ ./io/pfopen.cpp\ ./io/pfclose.cpp\ @@ -1240,4 +1229,6 @@ ./modules/Orthx/Orthx.h\ ./modules/Orthx/Orthx.cpp\ + ./modules/UpdateSpcsx/UpdateSpcsx.h\ + ./modules/UpdateSpcsx/UpdateSpcsx.cpp\ ./modules/Gradjx/Gradjx.h\ ./modules/Gradjx/Gradjx.cpp\ @@ -1270,7 +1261,4 @@ ./modules/ComputeStrainRatex/ComputeStrainRatex.h\ ./modules/ComputeStrainRatex/ComputeStrainRatex.cpp\ - ./modules/BuildNodeSetsx/BuildNodeSetsx.h\ - ./modules/BuildNodeSetsx/BuildNodeSetsx.cpp\ - ./modules/BuildNodeSetsx/PartitionSets.cpp\ ./modules/SpcNodesx/SpcNodesx.h\ ./modules/SpcNodesx/SpcNodesx.cpp\ @@ -1308,8 +1296,4 @@ ./modules/Reducevectorgtofx/Reducevectorgtofx.cpp\ ./modules/Reducevectorgtofx/Reducevectorgtofx.h\ - ./modules/Reducematrixfromgtofx/Reducematrixfromgtofx.cpp\ - ./modules/Reducematrixfromgtofx/Reducematrixfromgtofx.h\ - ./modules/Reduceloadfromgtofx/Reduceloadfromgtofx.h\ - ./modules/Reduceloadfromgtofx/Reduceloadfromgtofx.cpp\ ./modules/Reduceloadx/Reduceloadx.h\ ./modules/Reduceloadx/Reduceloadx.cpp\ @@ -1339,4 +1323,6 @@ ./modules/Solverx/Solverx.h\ ./modules/Solverx/DofTypesToIndexSet.cpp\ + ./modules/VecMergex/VecMergex.cpp\ + ./modules/VecMergex/VecMergex.h\ ./modules/Scotchx/Scotchx.cpp\ ./modules/Scotchx/Scotchx.h\ Index: /issm/trunk/src/c/modules/EnumToStringx/EnumToStringx.cpp =================================================================== --- /issm/trunk/src/c/modules/EnumToStringx/EnumToStringx.cpp (revision 8802) +++ /issm/trunk/src/c/modules/EnumToStringx/EnumToStringx.cpp (revision 8803) @@ -111,4 +111,6 @@ case PentaVertexInputEnum : return "PentaVertexInput"; case TriaVertexInputEnum : return "TriaVertexInput"; + case TriaVertexForcingEnum : return "TriaVertexForcing"; + case PentaVertexForcingEnum : return "PentaVertexForcing"; case ControlInputEnum : return "ControlInput"; case ParamEnum : return "Param"; Index: /issm/trunk/src/c/modules/StringToEnumx/StringToEnumx.cpp =================================================================== --- /issm/trunk/src/c/modules/StringToEnumx/StringToEnumx.cpp (revision 8802) +++ /issm/trunk/src/c/modules/StringToEnumx/StringToEnumx.cpp (revision 8803) @@ -109,4 +109,6 @@ else if (strcmp(name,"PentaVertexInput")==0) return PentaVertexInputEnum; else if (strcmp(name,"TriaVertexInput")==0) return TriaVertexInputEnum; + else if (strcmp(name,"TriaVertexForcing")==0) return TriaVertexForcingEnum; + else if (strcmp(name,"PentaVertexForcing")==0) return PentaVertexForcingEnum; else if (strcmp(name,"ControlInput")==0) return ControlInputEnum; else if (strcmp(name,"Param")==0) return ParamEnum; Index: /issm/trunk/src/c/solutions/ResetBoundaryConditions.cpp =================================================================== --- /issm/trunk/src/c/solutions/ResetBoundaryConditions.cpp (revision 8802) +++ /issm/trunk/src/c/solutions/ResetBoundaryConditions.cpp (revision 8803) @@ -9,8 +9,10 @@ void ResetBoundaryConditions(FemModel* femmodel, int analysis_type){ + /*variables: */ Vec yg=NULL; Vec ys=NULL; - int analysis_counter; - + Nodes* nodes=NULL; + int i; + _printf_(VerboseSolution(),"%s\n"," updating boundary conditions..."); @@ -18,17 +20,16 @@ femmodel->SetCurrentConfiguration(analysis_type); + /*recover nodes: */ + nodes=femmodel->nodes; + + /*retrieve boundary conditions from element inputs :*/ GetSolutionFromInputsx( &yg, femmodel->elements,femmodel->nodes, femmodel->vertices,femmodel->loads, femmodel->materials, femmodel->parameters); + Reducevectorgtosx(&ys,yg,femmodel->nodes,femmodel->parameters); - /*For this analysis_type, free existing boundary condition vector: */ - analysis_counter=femmodel->analysis_counter; - VecFree(&femmodel->m_ys[analysis_counter]); - - //Reduce from g to s set - Reducevectorgtosx(&ys,yg,femmodel->m_nodesets[analysis_counter],femmodel->parameters); - - /*Plug into femmodel->m_ys: */ - femmodel->m_ys[analysis_counter]=ys; + /*update spcs using this new vector of constraints: */ + UpdateSpcsx(femmodel->nodes,femmodel->parameters,ys); /*Free ressources:*/ + VecFree(&ys); VecFree(&yg); } Index: /issm/trunk/src/c/solutions/diagnostic_core.cpp =================================================================== --- /issm/trunk/src/c/solutions/diagnostic_core.cpp (revision 8802) +++ /issm/trunk/src/c/solutions/diagnostic_core.cpp (revision 8803) @@ -25,4 +25,5 @@ int solution_type; + /* recover parameters:*/ femmodel->parameters->FindParam(&dim,DimEnum); @@ -45,5 +46,5 @@ if(ishutter) surfaceslope_core(femmodel); if(isstokes) bedslope_core(femmodel); - + if(ishutter){ @@ -65,5 +66,5 @@ solver_nonlinear(femmodel,modify_loads); } - + if (ismacayealpattyn && isstokes){ Index: /issm/trunk/src/c/solvers/solver_adjoint_linear.cpp =================================================================== --- /issm/trunk/src/c/solvers/solver_adjoint_linear.cpp (revision 8802) +++ /issm/trunk/src/c/solvers/solver_adjoint_linear.cpp (revision 8803) @@ -13,26 +13,19 @@ /*intermediary: */ - Mat Kgg = NULL, Kff = NULL, Kfs = NULL; + Mat Kff = NULL, Kfs = NULL; Vec ug = NULL, uf = NULL; - Vec pg = NULL, pf = NULL; - Vec dg = NULL, df = NULL; - bool kffpartitioning = false; + Vec pf = NULL; + Vec df = NULL; + Vec ys = NULL; + int analysis_type; /*Recover parameters: */ - femmodel->parameters->FindParam(&kffpartitioning,KffEnum); + femmodel->parameters->FindParam(&analysis_type,AnalysisTypeEnum); - if(kffpartitioning){ - SystemMatricesx(NULL,&Kff, &Kfs, NULL,&pf, NULL, &df, NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); - Reduceloadx(pf, Kfs, femmodel->ys,true); MatFree(&Kfs); //true means spc = 0 - } - else{ - SystemMatricesx(&Kgg, NULL, NULL, &pg,NULL, &dg, NULL, NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); - Reducematrixfromgtofx(&Kff,&Kfs,Kgg,femmodel->nodesets); MatFree(&Kgg); - Reduceloadfromgtofx(&pf, pg, Kfs, femmodel->ys, femmodel->nodesets,true);VecFree(&pg); MatFree(&Kfs);//true means spc = 0 - Reducevectorgtofx(&df, dg, femmodel->nodesets,femmodel->parameters); VecFree(&dg); - } - + SystemMatricesx(&Kff, &Kfs, &pf, &df, NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); + CreateNodalConstraintsx(&ys,femmodel->nodes,analysis_type); + Reduceloadx(pf, Kfs, ys,true); MatFree(&Kfs); //true means spc = 0 Solverx(&uf, Kff, pf, NULL, df, femmodel->parameters); MatFree(&Kff); VecFree(&pf); VecFree(&df); - Mergesolutionfromftogx(&ug, uf,femmodel->ys,femmodel->nodesets,true); VecFree(&uf);//true means spc0 - InputUpdateFromSolutionx(femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,ug); VecFree(&ug); VecFree(&uf); + Mergesolutionfromftogx(&ug, uf,ys,femmodel->nodes,femmodel->parameters,true); VecFree(&uf);VecFree(&ys); //true means spc0 + InputUpdateFromSolutionx(femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,ug); VecFree(&ug); VecFree(&uf); } Index: /issm/trunk/src/c/solvers/solver_linear.cpp =================================================================== --- /issm/trunk/src/c/solvers/solver_linear.cpp (revision 8802) +++ /issm/trunk/src/c/solvers/solver_linear.cpp (revision 8803) @@ -11,26 +11,20 @@ /*intermediary: */ - Mat Kgg = NULL, Kff = NULL, Kfs = NULL; - Vec ug = NULL, uf = NULL; - Vec pg = NULL, pf = NULL; - Vec dg = NULL, df = NULL; - bool kffpartitioning; + Mat Kff = NULL, Kfs = NULL; + Vec ug = NULL; + Vec uf = NULL; + Vec pf = NULL; + Vec df = NULL; + Vec ys = NULL; + int analysis_type; /*Recover parameters: */ - femmodel->parameters->FindParam(&kffpartitioning,KffEnum); + femmodel->parameters->FindParam(&analysis_type,AnalysisTypeEnum); - if(kffpartitioning){ - SystemMatricesx(NULL,&Kff, &Kfs, NULL,&pf, NULL, &df, NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); - Reduceloadx(pf, Kfs, femmodel->ys); MatFree(&Kfs); - } - else{ - SystemMatricesx(&Kgg, NULL, NULL, &pg,NULL, &dg, NULL, NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); - Reducematrixfromgtofx(&Kff,&Kfs,Kgg,femmodel->nodesets); MatFree(&Kgg); - Reduceloadfromgtofx(&pf, pg, Kfs, femmodel->ys, femmodel->nodesets);VecFree(&pg); MatFree(&Kfs); - Reducevectorgtofx(&df, dg, femmodel->nodesets,femmodel->parameters); VecFree(&dg); - } - + SystemMatricesx(&Kff, &Kfs, &pf, &df, NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); + CreateNodalConstraintsx(&ys,femmodel->nodes,analysis_type); + Reduceloadx(pf, Kfs, ys); MatFree(&Kfs); Solverx(&uf, Kff, pf, NULL, df, femmodel->parameters); MatFree(&Kff); VecFree(&pf); VecFree(&df); - Mergesolutionfromftogx(&ug, uf,femmodel->ys,femmodel->nodesets);VecFree(&uf); + Mergesolutionfromftogx(&ug, uf,ys,femmodel->nodes,femmodel->parameters);VecFree(&uf);VecFree(&ys); InputUpdateFromSolutionx( femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,ug); VecFree(&ug); VecFree(&uf); } Index: /issm/trunk/src/c/solvers/solver_nonlinear.cpp =================================================================== --- /issm/trunk/src/c/solvers/solver_nonlinear.cpp (revision 8802) +++ /issm/trunk/src/c/solvers/solver_nonlinear.cpp (revision 8803) @@ -13,8 +13,10 @@ /*intermediary: */ - Mat Kgg = NULL, Kff = NULL, Kfs = NULL; - Vec ug = NULL, uf = NULL, old_ug= NULL, old_uf = NULL; - Vec pg = NULL, pf = NULL; - Vec dg = NULL, df = NULL; + Mat Kff = NULL, Kfs = NULL; + Vec ug = NULL, old_ug = NULL; + Vec uf = NULL, old_uf = NULL; + Vec pf = NULL; + Vec df = NULL; + Vec ys = NULL; Loads* loads=NULL; @@ -25,13 +27,14 @@ /*parameters:*/ - bool kffpartitioning=false; int min_mechanical_constraints; int max_nonlinear_iterations; + int analysis_type; + /*Recover parameters: */ - femmodel->parameters->FindParam(&kffpartitioning,KffEnum); femmodel->parameters->FindParam(&min_mechanical_constraints,MinMechanicalConstraintsEnum); femmodel->parameters->FindParam(&max_nonlinear_iterations,MaxNonlinearIterationsEnum); - + femmodel->parameters->FindParam(&analysis_type,AnalysisTypeEnum); + /*Were loads requested as output? : */ if(conserve_loads) loads=(Loads*)femmodel->loads->Copy(); //protect loads from being modified by the solution @@ -43,5 +46,5 @@ /*Start non-linear iteration using input velocity: */ GetSolutionFromInputsx(&ug, femmodel->elements, femmodel->nodes, femmodel->vertices, loads, femmodel->materials, femmodel->parameters); - Reducevectorgtofx(&uf, ug, femmodel->nodesets,femmodel->parameters); + Reducevectorgtofx(&uf, ug, femmodel->nodes,femmodel->parameters); //Update once again the solution to make sure that vx and vxold are similar (for next step in transient or steadystate) @@ -54,17 +57,9 @@ VecFree(&old_uf);old_uf=uf; - if(kffpartitioning){ - SystemMatricesx(NULL,&Kff, &Kfs, NULL,&pf, NULL, &df, NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,loads,femmodel->materials,femmodel->parameters); - Reduceloadx(pf, Kfs, femmodel->ys); MatFree(&Kfs); - } - else{ - SystemMatricesx(&Kgg, NULL, NULL, &pg,NULL, &dg, NULL, NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,loads,femmodel->materials,femmodel->parameters); - Reducematrixfromgtofx(&Kff,&Kfs,Kgg,femmodel->nodesets); MatFree(&Kgg); - Reduceloadfromgtofx(&pf, pg, Kfs, femmodel->ys, femmodel->nodesets); VecFree(&pg); MatFree(&Kfs); - Reducevectorgtofx(&df, dg, femmodel->nodesets,femmodel->parameters); VecFree(&dg); - } - + SystemMatricesx(&Kff, &Kfs, &pf, &df, NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,loads,femmodel->materials,femmodel->parameters); + CreateNodalConstraintsx(&ys,femmodel->nodes,analysis_type); + Reduceloadx(pf, Kfs, ys); MatFree(&Kfs); Solverx(&uf, Kff, pf, old_uf, df, femmodel->parameters); - Mergesolutionfromftogx(&ug, uf,femmodel->ys,femmodel->nodesets); + Mergesolutionfromftogx(&ug, uf,ys,femmodel->nodes,femmodel->parameters);VecFree(&ys); InputUpdateFromSolutionx( femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,ug); @@ -95,6 +90,4 @@ if(conserve_loads) delete loads; VecFree(&uf); - VecFree(&ug); VecFree(&old_uf); - VecFree(&old_ug); } Index: /issm/trunk/src/c/solvers/solver_stokescoupling_nonlinear.cpp =================================================================== --- /issm/trunk/src/c/solvers/solver_stokescoupling_nonlinear.cpp (revision 8802) +++ /issm/trunk/src/c/solvers/solver_stokescoupling_nonlinear.cpp (revision 8803) @@ -13,12 +13,13 @@ /*intermediary: */ - Mat Kgg_horiz = NULL, Kff_horiz = NULL, Kfs_horiz = NULL; - Vec ug_horiz = NULL, uf_horiz = NULL, old_ug_horiz= NULL, old_uf_horiz = NULL; - Vec pg_horiz = NULL, pf_horiz = NULL; - Vec dg_horiz = NULL, df_horiz = NULL; - Mat Kgg_vert = NULL, Kff_vert = NULL, Kfs_vert = NULL; + Mat Kff_horiz = NULL, Kfs_horiz = NULL; + Vec ug_horiz = NULL, uf_horiz = NULL, old_uf_horiz = NULL; + Vec pf_horiz = NULL; + Vec df_horiz = NULL; + Mat Kff_vert = NULL, Kfs_vert = NULL; Vec ug_vert = NULL, uf_vert = NULL; - Vec pg_vert = NULL, pf_vert = NULL; - Vec dg_vert = NULL, df_vert = NULL; + Vec pf_vert = NULL; + Vec df_vert = NULL; + Vec ys = NULL; int converged; int constraints_converged; @@ -27,10 +28,9 @@ /*parameters:*/ - bool kffpartitioning=false; int min_mechanical_constraints; int max_nonlinear_iterations; + int analysis_type; /*Recover parameters: */ - femmodel->parameters->FindParam(&kffpartitioning,KffEnum); femmodel->parameters->FindParam(&min_mechanical_constraints,MinMechanicalConstraintsEnum); femmodel->parameters->FindParam(&max_nonlinear_iterations,MaxNonlinearIterationsEnum); @@ -42,5 +42,5 @@ femmodel->SetCurrentConfiguration(DiagnosticHorizAnalysisEnum); GetSolutionFromInputsx(&ug_horiz, femmodel->elements, femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters); - Reducevectorgtofx(&uf_horiz, ug_horiz, femmodel->nodesets,femmodel->parameters); + Reducevectorgtofx(&uf_horiz, ug_horiz, femmodel->nodes,femmodel->parameters); for(;;){ @@ -49,25 +49,18 @@ /*First diagnostic horiz:*/ femmodel->SetCurrentConfiguration(DiagnosticHorizAnalysisEnum); - + femmodel->parameters->FindParam(&analysis_type,AnalysisTypeEnum); + //Update once again the solution to make sure that vx and vxold are similar (for next step in transient or steadystate) InputUpdateFromSolutionx(femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,ug_horiz); //save pointer to old velocity - VecFree(&old_ug_horiz);old_ug_horiz=ug_horiz; VecFree(&old_uf_horiz);old_uf_horiz=uf_horiz; - if(kffpartitioning){ - SystemMatricesx(NULL,&Kff_horiz, &Kfs_horiz, NULL,&pf_horiz, NULL, &df_horiz, NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); - Reduceloadx(pf_horiz, Kfs_horiz, femmodel->ys); MatFree(&Kfs_horiz); - } - else{ - SystemMatricesx(&Kgg_horiz, NULL, NULL, &pg_horiz,NULL, &dg_horiz,NULL,NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); - Reducematrixfromgtofx(&Kff_horiz,&Kfs_horiz,Kgg_horiz,femmodel->nodesets); MatFree(&Kgg_horiz); - Reduceloadfromgtofx(&pf_horiz, pg_horiz, Kfs_horiz, femmodel->ys, femmodel->nodesets); VecFree(&pg_horiz); MatFree(&Kfs_horiz); - Reducevectorgtofx(&df_horiz, dg_horiz, femmodel->nodesets,femmodel->parameters); VecFree(&dg_horiz); - } - + /*solve: */ + SystemMatricesx(&Kff_horiz, &Kfs_horiz, &pf_horiz, &df_horiz, NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); + CreateNodalConstraintsx(&ys,femmodel->nodes,analysis_type); + Reduceloadx(pf_horiz, Kfs_horiz, ys); MatFree(&Kfs_horiz); Solverx(&uf_horiz, Kff_horiz, pf_horiz, old_uf_horiz, df_horiz,femmodel->parameters); - Mergesolutionfromftogx(&ug_horiz, uf_horiz,femmodel->ys,femmodel->nodesets); + Mergesolutionfromftogx(&ug_horiz, uf_horiz,ys,femmodel->nodes,femmodel->parameters); VecFree(&ys); InputUpdateFromSolutionx( femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,ug_horiz); @@ -76,17 +69,12 @@ /*Second compute vertical velocity: */ femmodel->SetCurrentConfiguration(DiagnosticVertAnalysisEnum); - if(kffpartitioning){ - SystemMatricesx(NULL,&Kff_vert, &Kfs_vert, NULL,&pf_vert, NULL, &df_vert,NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); - Reduceloadx(pf_vert, Kfs_vert, femmodel->ys); MatFree(&Kfs_vert); - } - else{ - SystemMatricesx(&Kgg_vert, NULL, NULL, &pg_vert,NULL, &dg_horiz, NULL, NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); - Reducematrixfromgtofx(&Kff_vert,&Kfs_vert,Kgg_vert,femmodel->nodesets); MatFree(&Kgg_vert); - Reduceloadfromgtofx(&pf_vert, pg_vert, Kfs_vert, femmodel->ys, femmodel->nodesets);VecFree(&pg_vert); MatFree(&Kfs_vert); - Reducevectorgtofx(&df_vert, dg_vert, femmodel->nodesets,femmodel->parameters); VecFree(&dg_vert); - } - + femmodel->parameters->FindParam(&analysis_type,AnalysisTypeEnum); + + /*solve: */ + SystemMatricesx(&Kff_vert, &Kfs_vert, &pf_vert, &df_vert,NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); + CreateNodalConstraintsx(&ys,femmodel->nodes,analysis_type); + Reduceloadx(pf_vert, Kfs_vert, ys); MatFree(&Kfs_vert); Solverx(&uf_vert, Kff_vert, pf_vert, NULL, df_vert,femmodel->parameters); MatFree(&Kff_vert); VecFree(&pf_vert); VecFree(&df_vert); - Mergesolutionfromftogx(&ug_vert, uf_vert,femmodel->ys,femmodel->nodesets);VecFree(&uf_vert); + Mergesolutionfromftogx(&ug_vert, uf_vert,ys,femmodel->nodes,femmodel->parameters);VecFree(&uf_vert); VecFree(&ys); InputUpdateFromSolutionx( femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,ug_vert); VecFree(&ug_vert); VecFree(&uf_vert); @@ -102,6 +90,5 @@ /*clean-up*/ VecFree(&uf_horiz); - VecFree(&ug_horiz); VecFree(&old_uf_horiz); - VecFree(&old_ug_horiz); + VecFree(&ys); } Index: /issm/trunk/src/c/solvers/solver_thermal_nonlinear.cpp =================================================================== --- /issm/trunk/src/c/solvers/solver_thermal_nonlinear.cpp (revision 8802) +++ /issm/trunk/src/c/solvers/solver_thermal_nonlinear.cpp (revision 8803) @@ -14,13 +14,11 @@ Vec tf=NULL; Vec tf_old=NULL; + Vec ys=NULL; double melting_offset; /*intermediary: */ - Mat Kgg=NULL; Mat Kff=NULL; Mat Kfs=NULL; - Vec pg=NULL; Vec pf=NULL; - Vec dg=NULL; Vec df=NULL; @@ -35,11 +33,12 @@ int kflag,pflag; bool lowmem=0; - bool kffpartitioning; + int analysis_type; + /*Recover parameters: */ kflag=1; pflag=1; - femmodel->parameters->FindParam(&kffpartitioning,KffEnum); femmodel->parameters->FindParam(&lowmem,LowmemEnum); femmodel->parameters->FindParam(&min_thermal_constraints,MinThermalConstraintsEnum); + femmodel->parameters->FindParam(&analysis_type,AnalysisTypeEnum); count=1; @@ -51,21 +50,11 @@ for(;;){ - if(kffpartitioning){ - SystemMatricesx(NULL,&Kff, &Kfs, NULL,&pf,NULL, &df, &melting_offset,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); - Reduceloadx(pf, Kfs, femmodel->ys); MatFree(&Kfs); - } - else{ - SystemMatricesx(&Kgg,NULL, NULL, &pg,NULL, &dg, NULL, &melting_offset,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); - Reducematrixfromgtofx(&Kff,&Kfs,Kgg,femmodel->nodesets); MatFree(&Kgg); - Reduceloadfromgtofx(&pf, pg, Kfs, femmodel->ys, femmodel->nodesets);VecFree(&pg); MatFree(&Kfs); - Reducevectorgtofx(&df, dg, femmodel->nodesets,femmodel->parameters); VecFree(&dg); - } - - VecFree(&tf); + SystemMatricesx(&Kff, &Kfs, &pf,&df, &melting_offset,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); + CreateNodalConstraintsx(&ys,femmodel->nodes,analysis_type); + Reduceloadx(pf, Kfs, ys); MatFree(&Kfs); VecFree(&tf); Solverx(&tf, Kff, pf,tf_old, df, femmodel->parameters); VecFree(&tf_old); VecDuplicatePatch(&tf_old,tf); MatFree(&Kff);VecFree(&pf);VecFree(&tg); VecFree(&df); - - Mergesolutionfromftogx(&tg, tf,femmodel->ys,femmodel->nodesets); + Mergesolutionfromftogx(&tg, tf,ys,femmodel->nodes,femmodel->parameters); VecFree(&ys); InputUpdateFromSolutionx(femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,tg); @@ -90,3 +79,4 @@ VecFree(&tf); VecFree(&tf_old); + VecFree(&ys); } Index: /issm/trunk/src/m/classes/version/7.6/model.m =================================================================== --- /issm/trunk/src/m/classes/version/7.6/model.m (revision 8802) +++ /issm/trunk/src/m/classes/version/7.6/model.m (revision 8803) @@ -351,5 +351,4 @@ %solution parameters - kff=0; loadresults=0; @@ -782,7 +781,4 @@ md.petscoptions=addoptions(md.petscoptions,DiagnosticVertAnalysisEnum,mumpsoptions); - %solution speed-up - md.kff=1; - %hydrology: from Johnson's 2002 thesis, section 3.5.4 md.hydro_n=.02; Index: /issm/trunk/src/m/classes/version/7.7/model.m =================================================================== --- /issm/trunk/src/m/classes/version/7.7/model.m (revision 8802) +++ /issm/trunk/src/m/classes/version/7.7/model.m (revision 8803) @@ -351,5 +351,4 @@ %solution parameters - kff=0; loadresults=0; @@ -854,6 +853,4 @@ md.petscoptions=addoptions(md.petscoptions,DiagnosticVertAnalysisEnum,mumpsoptions); - %solution speed-up - md.kff=1; %hydrology: from Johnson's 2002 thesis, section 3.5.4 Index: /issm/trunk/src/m/model/marshall.m =================================================================== --- /issm/trunk/src/m/model/marshall.m (revision 8802) +++ /issm/trunk/src/m/model/marshall.m (revision 8803) @@ -193,5 +193,4 @@ WriteData(fid,md.shelf_dampening,'Integer','shelf_dampening'); WriteData(fid,md.waitonlock,'Scalar','waitonlock'); -WriteData(fid,md.kff,'Integer','kff'); %Thermal parameters Index: /issm/trunk/src/m/solutions/NewFemModel.m =================================================================== --- /issm/trunk/src/m/solutions/NewFemModel.m (revision 8802) +++ /issm/trunk/src/m/solutions/NewFemModel.m (revision 8803) @@ -21,6 +21,4 @@ %Initialize some fiels with empty celils - femmodel.m_nodesets=cell(nummodels,1); - femmodel.m_ys=cell(nummodels,1); femmodel.results=struct([]); @@ -44,10 +42,4 @@ [femmodel.nodes]=NodesDof(femmodel.nodes,femmodel.parameters); - issmprintf(VerboseMProcessor(),'%s',' create nodal constraints vector'); - femmodel.m_ys{i}=CreateNodalConstraints(femmodel.nodes,analysis_type); - - issmprintf(VerboseMProcessor(),'%s',' create node sets'); - femmodel.m_nodesets{i}=BuildNodeSets(femmodel.nodes,analysis_type); - issmprintf(VerboseMProcessor(),'%s',' configuring elements and loads'); [femmodel.elements,femmodel.loads,femmodel.nodes,femmodel.parameters] = ConfigureObjects( femmodel.elements, femmodel.loads, femmodel.nodes, femmodel.vertices,femmodel.materials,femmodel.parameters); Index: /issm/trunk/src/m/solutions/ResetBoundaryConditions.m =================================================================== --- /issm/trunk/src/m/solutions/ResetBoundaryConditions.m (revision 8802) +++ /issm/trunk/src/m/solutions/ResetBoundaryConditions.m (revision 8803) @@ -12,14 +12,8 @@ femmodel=SetCurrentConfiguration(femmodel,analysis_type); + %recreate constraints vector for the current analysis ug=GetSolutionFromInputs(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters); + ys=Reducevectorgtos(ug,femmode.nodes,femmodel.parameters); - %For this analysis_type, free existing boundary condition vectors: - analysis_counter=femmodel.parameters.AnalysisCounter+1; %matlab indexing on counter - femmodel.m_ys{analysis_counter}=[]; - - %Reduce from g to s set - ys=Reducevectorgtos(ug,femmodel.m_nodesets{analysis_counter},femmodel.parameters); - - %in the s-set - femmodel.m_ys{analysis_counter}=ys; - + %update values of spcs inside nodes + femmodel.nodes=UpdateSpcs(femmodel.nodes,femmodel.parameters,ys); Index: /issm/trunk/src/m/solutions/SetCurrentConfiguration.m =================================================================== --- /issm/trunk/src/m/solutions/SetCurrentConfiguration.m (revision 8802) +++ /issm/trunk/src/m/solutions/SetCurrentConfiguration.m (revision 8803) @@ -37,8 +37,4 @@ end - %activate matrices and vectors: - femmodel.nodesets=femmodel.m_nodesets{found}; - femmodel.ys=femmodel.m_ys{found}; - %Now, plug analysis_counter and analysis_type inside the parameters: %set counter and analyse_type Index: /issm/trunk/src/m/solvers/solver_adjoint_linear.m =================================================================== --- /issm/trunk/src/m/solvers/solver_adjoint_linear.m (revision 8802) +++ /issm/trunk/src/m/solvers/solver_adjoint_linear.m (revision 8803) @@ -6,20 +6,13 @@ %Get parameters - kffpartitioning=femmodel.parameters.Kff; + analysis_type=femmodel.parameters.AnalysisType; - if kffpartitioning, - [K_gg,K_ff,K_fs,p_g,p_f,dg,df,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters); - p_f = Reduceload( p_f, K_fs, femmodel.ys,true); - - else - [K_gg,K_ff,K_fs,p_g,p_f,dg,df,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters); - [K_ff, K_fs] = Reducematrixfromgtof( K_gg, femmodel.nodesets); - p_f = Reduceloadfromgtof( p_g, K_fs, femmodel.ys, femmodel.nodesets,true); - df=Reducevectorgtof( dg, femmodel.nodesets,femmodel.parameters); - end + [K_ff,K_fs,p_f,df,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters); + ys=CreateNodalConstraints(femmodel.nodes,analysis_type); + p_f = Reduceload( p_f, K_fs, ys,true); issmprintf(VerboseSolver(),'%s%g',' condition number of stiffness matrix: ',condest(K_ff)); u_f=Solver(K_ff,p_f,[],df,femmodel.parameters); - u_g= Mergesolutionfromftog( u_f, femmodel.ys, femmodel.nodesets,true); + u_g= Mergesolutionfromftog( u_f, femmodel.ys, femmodel.nodes,femmodel.parameters,true); [femmodel.elements,femmodel.materials]=InputUpdateFromSolution(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters,u_g); Index: /issm/trunk/src/m/solvers/solver_linear.m =================================================================== --- /issm/trunk/src/m/solvers/solver_linear.m (revision 8802) +++ /issm/trunk/src/m/solvers/solver_linear.m (revision 8803) @@ -6,19 +6,12 @@ %Get parameters - kffpartitioning=femmodel.parameters.Kff; + analysis_type=femmodel.parameters.AnalysisType; - if kffpartitioning, - [K_gg,K_ff,K_fs,p_g,p_f,dg,df,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters); - p_f = Reduceload( p_f, K_fs, femmodel.ys); - else - [K_gg,K_ff,K_fs,p_g,p_f,dg,df,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters); - [K_ff, K_fs] = Reducematrixfromgtof( K_gg, femmodel.nodesets); - p_f = Reduceloadfromgtof( p_g, K_fs, femmodel.ys, femmodel.nodesets); - df=Reducevectorgtof( dg, femmodel.nodesets,femmodel.parameters); - end + [K_ff,K_fs,p_f,df,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters); + ys=CreateNodalConstraints(femmodel.nodes,analysis_type); + p_f = Reduceload( p_f, K_fs, ys); issmprintf(VerboseSolver(),'%s%g',' condition number of stiffness matrix: ',condest(K_ff)); u_f=Solver(K_ff,p_f,[],df,femmodel.parameters); - u_g= Mergesolutionfromftog( u_f, femmodel.ys, femmodel.nodesets); - + u_g= Mergesolutionfromftog( u_f, ys, femmodel.nodes,femmodel.parameters); [femmodel.elements,femmodel.materials]=InputUpdateFromSolution(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters,u_g); Index: /issm/trunk/src/m/solvers/solver_nonlinear.m =================================================================== --- /issm/trunk/src/m/solvers/solver_nonlinear.m (revision 8802) +++ /issm/trunk/src/m/solvers/solver_nonlinear.m (revision 8803) @@ -6,6 +6,6 @@ %Branch on partitioning schema requested - kffpartitioning=femmodel.parameters.Kff; min_mechanical_constraints=femmodel.parameters.MinMechanicalConstraints; + analysis_type=femmodel.parameters.AnalysisType; %keep a copy of loads for now @@ -17,5 +17,5 @@ %Start non-linear iteration using input velocity: ug=GetSolutionFromInputs(femmodel.elements, femmodel.nodes, femmodel.vertices, loads, femmodel.materials, femmodel.parameters); - uf=Reducevectorgtof( ug, femmodel.nodesets,femmodel.parameters); + uf=Reducevectorgtof( ug, femmodel.nodes,femmodel.parameters); %Update the solution to make sure that vx and vxold are similar @@ -28,16 +28,10 @@ old_uf=uf; - if kffpartitioning, - [K_gg,K_ff,K_fs,p_g,p_f,dg,df,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,loads,femmodel.materials,femmodel.parameters); - p_f = Reduceload( p_f, K_fs, femmodel.ys); - else - [K_gg,K_ff,K_fs,p_g,p_f,dg,df,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,loads,femmodel.materials,femmodel.parameters); - [K_ff, K_fs] = Reducematrixfromgtof( K_gg, femmodel.nodesets); - p_f = Reduceloadfromgtof( p_g, K_fs, femmodel.ys, femmodel.nodesets); - df=Reducevectorgtof( dg, femmodel.nodesets,femmodel.parameters); - end + [K_ff,K_fs,p_f,df,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,loads,femmodel.materials,femmodel.parameters); + ys=CreateNodalConstraints(femmodel.nodes,analysis_type); + p_f = Reduceload( p_f, K_fs, femmodel.ys); uf=Solver(K_ff,p_f,old_uf,df,femmodel.parameters); - ug= Mergesolutionfromftog( uf, femmodel.ys, femmodel.nodesets); + ug= Mergesolutionfromftog( uf, ys, femmodel.nodes,femmodel.parameters); [femmodel.elements,femmodel.materials]=InputUpdateFromSolution(femmodel.elements,femmodel.nodes,femmodel.vertices,loads,femmodel.materials,femmodel.parameters,ug); Index: /issm/trunk/src/m/solvers/solver_stokescoupling_nonlinear.m =================================================================== --- /issm/trunk/src/m/solvers/solver_stokescoupling_nonlinear.m (revision 8802) +++ /issm/trunk/src/m/solvers/solver_stokescoupling_nonlinear.m (revision 8803) @@ -4,7 +4,4 @@ % Usage: % [femmodel]=solver_couplingstokes_nonlinear(femmodel,conserve_loads) - - %Get parameters - kffpartitioning=femmodel.parameters.Kff; %initialize solution vector @@ -14,5 +11,5 @@ femmodel=SetCurrentConfiguration(femmodel,DiagnosticHorizAnalysisEnum); ug_horiz=GetSolutionFromInputs(femmodel.elements, femmodel.nodes, femmodel.vertices, femmodel.loads, femmodel.materials, femmodel.parameters); - uf_horiz=Reducevectorgtof( ug_horiz, femmodel.nodesets,femmodel.parameters); + uf_horiz=Reducevectorgtof( ug_horiz, femmodel.nodes,femmodel.parameters); while(~converged), @@ -20,4 +17,5 @@ %First compute the horizontal velocity femmodel=SetCurrentConfiguration(femmodel,DiagnosticHorizAnalysisEnum); + analysis_type=femmodel.parameters.AnalysisType; %Update the solution to make sure that vx and vxold are similar @@ -28,17 +26,10 @@ old_uf_horiz=uf_horiz; - if kffpartitioning, - [K_gg_horiz,K_ff_horiz,K_fs_horiz,p_g_horiz,p_f_horiz,d_g_horiz,d_f_horiz,kmax_horiz]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters); - p_f_horiz = Reduceload( p_f_horiz, K_fs_horiz, femmodel.ys); - else - [K_gg_horiz,K_ff_horiz,K_fs_horiz,p_g_horiz,p_f_horiz,d_g_horiz,d_f_horiz,kmax]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters); - [K_ff_horiz, K_fs_horiz] = Reducematrixfromgtof( K_gg_horiz, femmodel.nodesets); - p_f_horiz = Reduceloadfromgtof( p_g_horiz, K_fs_horiz, femmodel.ys, femmodel.nodesets); - d_f_horiz=Reducevectorgtof( d_g_horiz, femmodel.nodesets,femmodel.parameters); - - end + [K_ff_horiz,K_fs_horiz,p_f_horiz,d_f_horiz,kmax_horiz]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters); + ys=CreateNodalConstraints(femmodel.nodes,analysis_type); + p_f_horiz = Reduceload( p_f_horiz, K_fs_horiz, ys); uf_horiz=Solver(K_ff_horiz,p_f_horiz,old_uf_horiz,d_f_horiz,femmodel.parameters); - ug_horiz= Mergesolutionfromftog( uf_horiz, femmodel.ys, femmodel.nodesets); + ug_horiz= Mergesolutionfromftog( uf_horiz, ys, femmodel.nodes,femmodel.parameters); [femmodel.elements,femmodel.materials]=InputUpdateFromSolution(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters,ug_horiz); @@ -49,17 +40,12 @@ %Then compute vertical velocity femmodel=SetCurrentConfiguration(femmodel,DiagnosticVertAnalysisEnum); + analysis_type=femmodel.parameters.AnalysisType; - if kffpartitioning, - [K_gg_vert,K_ff_vert,K_fs_vert,p_g_vert,p_f_vert,d_g_vert,d_f_vert,kmax_vert]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters); - p_f_vert = Reduceload( p_f_vert, K_fs_vert, femmodel.ys); - else - [K_gg_vert,K_ff_vert,K_fs_vert,p_g_vert,p_f_vert,d_g_vert,d_f_vert,kmax_vert]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters); - [K_ff_vert, K_fs_vert] = Reducematrixfromgtof( K_gg_vert, femmodel.nodesets); - p_f_vert = Reduceloadfromgtof( p_g_vert, K_fs_vert, femmodel.ys, femmodel.nodesets); - d_f_vert=Reducevectorgtof( d_g_vert, femmodel.nodesets,femmodel.parameters); - end + [K_ff_vert,K_fs_vert,p_f_vert,d_f_vert,kmax_vert]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters); + ys=CreateNodalConstraints(femmodel.nodes,analysis_type); + p_f_vert = Reduceload( p_f_vert, K_fs_vert, ys); uf_vert=Solver(K_ff_vert,p_f_vert,[],d_f_vert,femmodel.parameters); - ug_vert= Mergesolutionfromftog( uf_vert, femmodel.ys, femmodel.nodesets); + ug_vert= Mergesolutionfromftog( uf_vert, ys, femmodel.nodes,femmodel.parameters); [femmodel.elements,femmodel.materials]=InputUpdateFromSolution(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters,ug_vert); Index: /issm/trunk/src/m/solvers/solver_thermal_nonlinear.m =================================================================== --- /issm/trunk/src/m/solvers/solver_thermal_nonlinear.m (revision 8802) +++ /issm/trunk/src/m/solvers/solver_thermal_nonlinear.m (revision 8803) @@ -10,5 +10,5 @@ %Get parameters - kffpartitioning=femmodel.parameters.Kff; + analysis_type=femmodel.parameters.AnalysisType; %stiffness and load generation only: @@ -23,18 +23,11 @@ while(~converged), - if kffpartitioning, - [K_gg,K_ff,K_fs,p_g,p_f,d_g,d_f,melting_offset]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters); - p_f = Reduceload( p_f, K_fs, femmodel.ys); - else - [K_gg,K_ff,K_fs,p_g,p_f,d_g,d_f,melting_offset]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters); - [K_ff, K_fs] = Reducematrixfromgtof( K_gg, femmodel.nodesets); - p_f = Reduceloadfromgtof( p_g, K_fs, femmodel.ys, femmodel.nodesets); - d_f=Reducevectorgtof( d_g, femmodel.nodesets,femmodel.parameters); - - end + [K_ff,K_fs,p_f,d_f,melting_offset]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters); + ys=CreateNodalConstraints(femmodel.nodes,analysis_type); + p_f = Reduceload( p_f, K_fs, ys); issmprintf(VerboseSolver(),'%s%g',' condition number of stiffness matrix: ',condest(K_ff)); t_f=Solver(K_ff,p_f,[],d_f,femmodel.parameters); - t_g= Mergesolutionfromftog( t_f, femmodel.ys, femmodel.nodesets); + t_g= Mergesolutionfromftog( t_f, ys, femmodel.nodes,femmodel.parameters); [femmodel.elements,femmodel.materials]=InputUpdateFromSolution(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters,t_g);