Changeset 5686
- Timestamp:
- 09/07/10 10:59:58 (15 years ago)
- Location:
- issm/trunk/src
- Files:
-
- 7 edited
Legend:
- Unmodified
- Added
- Removed
-
issm/trunk/src/c/modules/Solverx/Solverx.cpp
r4954 r5686 4 4 5 5 #include "./Solverx.h" 6 7 6 #include "../../shared/shared.h" 8 7 #include "../../include/include.h" … … 14 13 #endif 15 14 16 17 void Solverx( Vec* puf, Mat Kff, Vec pf, Vec uf0, char* solver_string){ 15 void Solverx(Vec* puf, Mat Kff, Vec pf, Vec uf0,Parameters* parameters){ 18 16 19 17 /*output: */ … … 29 27 PetscTruth flag; 30 28 int solver_type; 29 char* solver_string=NULL; 31 30 32 31 /*First, check that f-set is not NULL, ie model is fully constrained: */ 33 32 if(!Kff){ 34 *puf=NULL; 35 return; 33 *puf=NULL; return; 36 34 } 37 35 … … 46 44 47 45 /*Before preparing the solver, add options to the options database*/ 46 parameters->FindParam(&solver_string,SolverStringEnum); 48 47 PetscOptionsInsertMultipleString(solver_string); 49 48 … … 95 94 } 96 95 } 97 98 96 KSPSolve(ksp,pf,uf); 99 97 100 98 /*Check convergence*/ 101 99 KSPGetIterationNumber(ksp,&iteration_number); 102 if (iteration_number<0){ 103 ISSMERROR("%s%i"," Solver diverged at iteration number: ",-iteration_number); 104 } 100 if (iteration_number<0) ISSMERROR("%s%i"," Solver diverged at iteration number: ",-iteration_number); 105 101 106 102 /*Free ressources:*/ 107 103 KSPFree(&ksp); 108 104 xfree((void**)&solver_string); 109 105 110 106 /*Assign output pointers:*/ -
issm/trunk/src/c/modules/Solverx/Solverx.h
r3913 r5686 8 8 #include "../../objects/objects.h" 9 9 10 11 10 /* local prototypes: */ 12 void Solverx( Vec* puf, Mat Kff, Vec pf, Vec uf0, char* solver_string);11 void Solverx( Vec* puf, Mat Kff, Vec pf, Vec uf0,Parameters* parameters); 13 12 14 13 #endif /* _SOLVERX_H */ -
issm/trunk/src/c/solvers/solver_adjoint_linear.cpp
r5057 r5686 8 8 #include "../modules/modules.h" 9 9 10 void solver_adjoint_linear(FemModel* fem ){10 void solver_adjoint_linear(FemModel* femmodel){ 11 11 /*This is axactly the same solver as solver_linear except that Reduceloadfromgtofx and Mergesolutionfromftogx 12 12 * use the flag "true" so that all spc are taken as 0*/ … … 15 15 int kflag,pflag; 16 16 int verbose=0; 17 char* solver_string=NULL;18 17 19 18 /*output: */ … … 30 29 /*Recover parameters: */ 31 30 kflag=1; pflag=1; 32 fem->parameters->FindParam(&verbose,VerboseEnum); 33 fem->parameters->FindParam(&solver_string,SolverStringEnum); 31 femmodel->parameters->FindParam(&verbose,VerboseEnum); 34 32 35 33 //*Generate system matrices 36 34 if(verbose)_printf_(" Generating matrices\n"); 37 SystemMatricesx(&Kgg, &pg,fem ->elements,fem->nodes,fem->vertices,fem->loads,fem->materials,fem->parameters,kflag,pflag);35 SystemMatricesx(&Kgg, &pg,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,kflag,pflag); 38 36 39 37 if(verbose)_printf_(" Generating penalty matrices\n"); 40 38 //*Generate penalty system matrices 41 PenaltySystemMatricesx(Kgg, pg,NULL,fem ->elements,fem->nodes,fem->vertices,fem->loads,fem->materials,fem->parameters,kflag,pflag);39 PenaltySystemMatricesx(Kgg, pg,NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,kflag,pflag); 42 40 43 41 /*!Reduce matrix from g to f size:*/ 44 42 if(verbose)_printf_(" reducing matrix from g to f set\n"); 45 Reducematrixfromgtofx(&Kff,&Kfs,Kgg,fem ->nodesets); MatFree(&Kgg);43 Reducematrixfromgtofx(&Kff,&Kfs,Kgg,femmodel->nodesets); MatFree(&Kgg); 46 44 47 45 /*!Reduce load from g to f size: */ 48 46 if(verbose)_printf_(" reducing load from g to f set\n"); //true means spc = 0 49 Reduceloadfromgtofx(&pf, pg, Kfs, fem ->ys, fem->nodesets,true);VecFree(&pg); MatFree(&Kfs);47 Reduceloadfromgtofx(&pf, pg, Kfs, femmodel->ys, femmodel->nodesets,true);VecFree(&pg); MatFree(&Kfs); 50 48 51 49 /*Solve: */ 52 50 if(verbose)_printf_(" solving\n"); 53 Solverx(&uf, Kff, pf, NULL, solver_string); MatFree(&Kff); VecFree(&pf);51 Solverx(&uf, Kff, pf, NULL, femmodel->parameters); MatFree(&Kff); VecFree(&pf); 54 52 55 53 //Merge back to g set 56 54 if(verbose)_printf_(" merging solution from f to g set\n");//true means spc=0 57 Mergesolutionfromftogx(&ug, uf,fem ->ys,fem->nodesets,true);VecFree(&uf);55 Mergesolutionfromftogx(&ug, uf,femmodel->ys,femmodel->nodesets,true);VecFree(&uf); 58 56 59 57 //Update inputs using new solution: 60 InputUpdateFromSolutionx( fem ->elements,fem->nodes, fem->vertices, fem->loads, fem->materials, fem->parameters,ug);58 InputUpdateFromSolutionx( femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,ug); 61 59 62 60 /*free ressources: */ 63 xfree((void**)&solver_string);64 61 VecFree(&ug); 65 62 VecFree(&uf); -
issm/trunk/src/c/solvers/solver_diagnostic_nonlinear.cpp
r5057 r5686 35 35 /*parameters:*/ 36 36 int kflag,pflag; 37 char* solver_string=NULL;38 37 int verbose=0; 39 38 int dim; … … 41 40 /*Recover parameters: */ 42 41 kflag=1; pflag=1; 43 femmodel->parameters->FindParam(&solver_string,SolverStringEnum);44 42 femmodel->parameters->FindParam(&dim,DimEnum); 45 43 femmodel->parameters->FindParam(&verbose,VerboseEnum); … … 91 89 92 90 if(verbose)_printf_(" solving\n"); 93 Solverx(&uf, Kff, pf, old_uf, solver_string);91 Solverx(&uf, Kff, pf, old_uf, femmodel->parameters); 94 92 95 93 if(verbose)_printf_(" merging solution from f to g set\n"); … … 136 134 VecFree(&old_uf); 137 135 VecFree(&old_ug); 138 xfree((void**)&solver_string);139 136 140 137 } -
issm/trunk/src/c/solvers/solver_linear.cpp
r5057 r5686 8 8 #include "../modules/modules.h" 9 9 10 void solver_linear(FemModel* fem ){10 void solver_linear(FemModel* femmodel){ 11 11 12 12 /*parameters:*/ 13 13 int kflag,pflag; 14 14 int verbose=0; 15 char* solver_string=NULL;16 15 17 16 /*output: */ … … 28 27 /*Recover parameters: */ 29 28 kflag=1; pflag=1; 30 fem->parameters->FindParam(&verbose,VerboseEnum); 31 fem->parameters->FindParam(&solver_string,SolverStringEnum); 29 femmodel->parameters->FindParam(&verbose,VerboseEnum); 32 30 33 31 //*Generate system matrices 34 32 if(verbose)_printf_(" Generating matrices\n"); 35 SystemMatricesx(&Kgg, &pg,fem ->elements,fem->nodes,fem->vertices,fem->loads,fem->materials,fem->parameters,kflag,pflag);33 SystemMatricesx(&Kgg, &pg,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,kflag,pflag); 36 34 37 35 if(verbose)_printf_(" Generating penalty matrices\n"); 38 36 //*Generate penalty system matrices 39 PenaltySystemMatricesx(Kgg, pg,NULL,fem ->elements,fem->nodes,fem->vertices,fem->loads,fem->materials,fem->parameters,kflag,pflag);37 PenaltySystemMatricesx(Kgg, pg,NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,kflag,pflag); 40 38 41 39 /*!Reduce matrix from g to f size:*/ 42 40 if(verbose)_printf_(" reducing matrix from g to f set\n"); 43 Reducematrixfromgtofx(&Kff,&Kfs,Kgg,fem ->nodesets); MatFree(&Kgg);41 Reducematrixfromgtofx(&Kff,&Kfs,Kgg,femmodel->nodesets); MatFree(&Kgg); 44 42 45 43 /*!Reduce load from g to f size: */ 46 44 if(verbose)_printf_(" reducing load from g to f set\n"); 47 Reduceloadfromgtofx(&pf, pg, Kfs, fem ->ys, fem->nodesets);VecFree(&pg); MatFree(&Kfs);45 Reduceloadfromgtofx(&pf, pg, Kfs, femmodel->ys, femmodel->nodesets);VecFree(&pg); MatFree(&Kfs); 48 46 49 47 /*Solve: */ 50 48 if(verbose)_printf_(" solving\n"); 51 Solverx(&uf, Kff, pf, NULL, solver_string); MatFree(&Kff); VecFree(&pf);49 Solverx(&uf, Kff, pf, NULL, femmodel->parameters); MatFree(&Kff); VecFree(&pf); 52 50 53 51 //Merge back to g set 54 52 if(verbose)_printf_(" merging solution from f to g set\n"); 55 Mergesolutionfromftogx(&ug, uf,fem ->ys,fem->nodesets);VecFree(&uf);53 Mergesolutionfromftogx(&ug, uf,femmodel->ys,femmodel->nodesets);VecFree(&uf); 56 54 57 55 //Update inputs using new solution: 58 InputUpdateFromSolutionx( fem ->elements,fem->nodes, fem->vertices, fem->loads, fem->materials, fem->parameters,ug);56 InputUpdateFromSolutionx( femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,ug); 59 57 60 58 /*free ressources: */ 61 xfree((void**)&solver_string);62 59 VecFree(&ug); 63 60 VecFree(&uf); -
issm/trunk/src/c/solvers/solver_thermal_nonlinear.cpp
r5057 r5686 8 8 #include "../modules/modules.h" 9 9 10 void solver_thermal_nonlinear(FemModel* fem ){10 void solver_thermal_nonlinear(FemModel* femmodel){ 11 11 12 12 /*solution : */ … … 34 34 /*parameters:*/ 35 35 int kflag,pflag; 36 char* solver_string=NULL;37 36 int verbose=0; 38 37 bool lowmem=0; … … 41 40 kflag=1; pflag=1; 42 41 43 fem->parameters->FindParam(&solver_string,SolverStringEnum); 44 fem->parameters->FindParam(&verbose,VerboseEnum); 45 fem->parameters->FindParam(&lowmem,LowmemEnum); 46 fem->parameters->FindParam(&min_thermal_constraints,MinThermalConstraintsEnum); 42 femmodel->parameters->FindParam(&verbose,VerboseEnum); 43 femmodel->parameters->FindParam(&lowmem,LowmemEnum); 44 femmodel->parameters->FindParam(&min_thermal_constraints,MinThermalConstraintsEnum); 47 45 48 46 count=1; … … 54 52 55 53 if(count==1) reset_penalties=1; else reset_penalties=0; 56 InputUpdateFromConstantx( fem ->elements,fem->nodes, fem->vertices, fem->loads, fem->materials, fem->parameters,reset_penalties,ResetPenaltiesEnum);54 InputUpdateFromConstantx( femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,reset_penalties,ResetPenaltiesEnum); 57 55 58 56 //*Generate system matrices … … 61 59 /*Compute Kgg_nopenalty and pg_nopenalty once for all: */ 62 60 if (count==1){ 63 SystemMatricesx(&Kgg_nopenalty, &pg_nopenalty,fem ->elements,fem->nodes,fem->vertices,fem->loads,fem->materials,fem->parameters,kflag,pflag);61 SystemMatricesx(&Kgg_nopenalty, &pg_nopenalty,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,kflag,pflag); 64 62 } 65 63 … … 69 67 70 68 //apply penalties each time 71 PenaltySystemMatricesx(Kgg, pg,&melting_offset,fem ->elements,fem->nodes,fem->vertices,fem->loads,fem->materials,fem->parameters,kflag,pflag);69 PenaltySystemMatricesx(Kgg, pg,&melting_offset,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,kflag,pflag); 72 70 } 73 71 else{ 74 SystemMatricesx(&Kgg, &pg,fem ->elements,fem->nodes,fem->vertices,fem->loads,fem->materials,fem->parameters,kflag,pflag);72 SystemMatricesx(&Kgg, &pg,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,kflag,pflag); 75 73 //apply penalties 76 PenaltySystemMatricesx(Kgg, pg,&melting_offset,fem ->elements,fem->nodes,fem->vertices,fem->loads,fem->materials,fem->parameters,kflag,pflag);74 PenaltySystemMatricesx(Kgg, pg,&melting_offset,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,kflag,pflag); 77 75 } 78 76 79 77 /*!Reduce matrix from g to f size:*/ 80 Reducematrixfromgtofx(&Kff,&Kfs,Kgg,fem ->nodesets);78 Reducematrixfromgtofx(&Kff,&Kfs,Kgg,femmodel->nodesets); 81 79 82 80 /*Free ressources: */ … … 84 82 85 83 if(verbose)_printf_(" reducing load from g to f set\n"); 86 Reduceloadfromgtofx(&pf, pg, Kfs, fem ->ys, fem->nodesets);84 Reduceloadfromgtofx(&pf, pg, Kfs, femmodel->ys, femmodel->nodesets); 87 85 88 86 //no need for pg and Kfs anymore … … 93 91 if(verbose)_printf_("%s\n","solving"); 94 92 VecFree(&tf); 95 Solverx(&tf, Kff, pf,tf_old, solver_string);93 Solverx(&tf, Kff, pf,tf_old, femmodel->parameters); 96 94 VecFree(&tf_old); VecDuplicatePatch(&tf_old,tf); 97 95 … … 100 98 101 99 if(verbose)_printf_(" merging solution from f to g set\n"); 102 Mergesolutionfromftogx(&tg, tf,fem ->ys,fem->nodesets);100 Mergesolutionfromftogx(&tg, tf,femmodel->ys,femmodel->nodesets); 103 101 104 102 //Update inputs using new solution: 105 103 if (verbose) _printf_(" updating inputs\n"); 106 InputUpdateFromSolutionx(fem ->elements,fem->nodes, fem->vertices, fem->loads, fem->materials, fem->parameters,tg);104 InputUpdateFromSolutionx(femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,tg); 107 105 108 106 //Deal with penalty loads 109 107 if(verbose)_printf_(" penalty constraints\n"); 110 PenaltyConstraintsx(&constraints_converged, &num_unstable_constraints, fem ->elements,fem->nodes,fem->vertices,fem->loads,fem->materials,fem->parameters);108 PenaltyConstraintsx(&constraints_converged, &num_unstable_constraints, femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); 111 109 112 110 if (!converged){ … … 120 118 121 119 //add melting_offset to inputs: 122 InputUpdateFromConstantx( fem ->elements,fem->nodes, fem->vertices, fem->loads, fem->materials, fem->parameters,melting_offset,MeltingOffsetEnum);120 InputUpdateFromConstantx( femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,melting_offset,MeltingOffsetEnum); 123 121 124 122 /*Free ressources: */ … … 128 126 VecFree(&tf); 129 127 VecFree(&tf_old); 130 xfree((void**)&solver_string);131 128 } -
issm/trunk/src/mex/Solver/Solver.cpp
r4526 r5686 37 37 FetchData(&pf,PF); 38 38 FetchData(&uf0,UF0); 39 }40 41 /*Run solver code: */42 if (strcmp(solver_string,matlabstring)!=0){43 39 /*Petsc solver: */ 44 Solverx(&uf, Kff, pf, uf0, solver_string); 40 Solverx(&uf, Kff, pf, uf0, parameters); 41 /*Write output*/ 42 WriteData(UF,uf); 45 43 } 46 44 else{ … … 51 49 } 52 50 53 /*write output datasets if running Petsc solver: */54 if (strcmp(solver_string,matlabstring)!=0) WriteData(UF,uf);55 56 51 /*Free ressources: */ 57 if (strcmp(solver_string,matlabstring)!=0){ 58 MatFree(&Kff); 59 VecFree(&pf); 60 VecFree(&uf0); 61 VecFree(&uf); 62 delete parameters; 63 } 52 MatFree(&Kff); 53 VecFree(&pf); 54 VecFree(&uf0); 55 VecFree(&uf); 56 delete parameters; 64 57 xfree((void**)&solver_string); 65 58 … … 71 64 { 72 65 _printf_("\n"); 73 _printf_(" usage: [uf] = %s(Kff,pf,uf0, solver_string);\n",__FUNCT__);66 _printf_(" usage: [uf] = %s(Kff,pf,uf0,parameters);\n",__FUNCT__); 74 67 _printf_("\n"); 75 68 }
Note:
See TracChangeset
for help on using the changeset viewer.