Changeset 6001
- Timestamp:
- 09/23/10 15:24:53 (14 years ago)
- Location:
- issm/trunk/src/c/solvers
- Files:
-
- 4 edited
Legend:
- Unmodified
- Added
- Removed
-
issm/trunk/src/c/solvers/solver_adjoint_linear.cpp
r5772 r6001 13 13 14 14 /*intermediary: */ 15 Mat Kgg = NULL, Kff=NULL, Kfs=NULL; 16 Vec ug = NULL, uf=NULL; 17 Vec pg = NULL, pf=NULL; 15 Mat Kgg = NULL, Kff = NULL, Kfs = NULL; 16 Vec ug = NULL, uf = NULL; 17 Vec pg = NULL, pf = NULL; 18 bool kffpartitioning = false; 18 19 19 SystemMatricesx(&Kgg, NULL, NULL, &pg,NULL, NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); 20 /*Recover parameters: */ 21 femmodel->parameters->FindParam(&kffpartitioning,KffEnum); 22 23 if(kffpartitioning){ 24 SystemMatricesx(NULL,&Kff, &Kfs, NULL,&pf, NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); 25 Reduceloadx(pf, Kfs, femmodel->ys,true); MatFree(&Kfs); //true means spc = 0 26 } 27 else{ 28 SystemMatricesx(&Kgg, NULL, NULL, &pg,NULL, NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); 29 Reducematrixfromgtofx(&Kff,&Kfs,Kgg,femmodel->nodesets,femmodel->parameters); MatFree(&Kgg); 30 Reduceloadfromgtofx(&pf, pg, Kfs, femmodel->ys, femmodel->nodesets,femmodel->parameters,true);VecFree(&pg); MatFree(&Kfs);//true means spc = 0 31 } 20 32 21 33 Reducematrixfromgtofx(&Kff,&Kfs,Kgg,femmodel->nodesets,femmodel->parameters); MatFree(&Kgg); 22 Reduceloadfromgtofx(&pf, pg, Kfs, femmodel->ys, femmodel->nodesets,femmodel->parameters,true);VecFree(&pg); MatFree(&Kfs);//true means spc = 023 34 24 35 Solverx(&uf, Kff, pf, NULL, femmodel->parameters); MatFree(&Kff); VecFree(&pf); -
issm/trunk/src/c/solvers/solver_diagnostic_nonlinear.cpp
r5997 r6001 24 24 /*parameters:*/ 25 25 int verbose=0; 26 bool kff =false;26 bool kffpartitioning=false; 27 27 int min_mechanical_constraints; 28 28 int max_nonlinear_iterations; … … 30 30 /*Recover parameters: */ 31 31 femmodel->parameters->FindParam(&verbose,VerboseEnum); 32 femmodel->parameters->FindParam(&kff ,KffEnum);32 femmodel->parameters->FindParam(&kffpartitioning,KffEnum); 33 33 femmodel->parameters->FindParam(&min_mechanical_constraints,MinMechanicalConstraintsEnum); 34 34 femmodel->parameters->FindParam(&max_nonlinear_iterations,MaxNonlinearIterationsEnum); … … 54 54 VecFree(&old_uf);old_uf=uf; 55 55 56 if(kff ){56 if(kffpartitioning){ 57 57 SystemMatricesx(NULL,&Kff, &Kfs, NULL,&pf, NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,loads,femmodel->materials,femmodel->parameters); 58 Reduceloadx(pf, Kfs, femmodel->ys , femmodel->parameters); MatFree(&Kfs);58 Reduceloadx(pf, Kfs, femmodel->ys); MatFree(&Kfs); 59 59 } 60 60 else{ -
issm/trunk/src/c/solvers/solver_linear.cpp
r5772 r6001 11 11 12 12 /*intermediary: */ 13 Mat Kgg = NULL, Kff = NULL, Kfs = NULL; 14 Vec ug = NULL, uf = NULL; 15 Vec pg = NULL, pf = NULL; 13 Mat Kgg = NULL, Kff = NULL, Kfs = NULL; 14 Vec ug = NULL, uf = NULL; 15 Vec pg = NULL, pf = NULL; 16 bool kffpartitioning; 16 17 17 SystemMatricesx(&Kgg,NULL, NULL, &pg,NULL, NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); 18 /*Recover parameters: */ 19 femmodel->parameters->FindParam(&kffpartitioning,KffEnum); 18 20 19 Reducematrixfromgtofx(&Kff,&Kfs,Kgg,femmodel->nodesets,femmodel->parameters); MatFree(&Kgg); 20 Reduceloadfromgtofx(&pf, pg, Kfs, femmodel->ys, femmodel->nodesets,femmodel->parameters);VecFree(&pg); MatFree(&Kfs); 21 if(kffpartitioning){ 22 SystemMatricesx(NULL,&Kff, &Kfs, NULL,&pf, NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); 23 Reduceloadx(pf, Kfs, femmodel->ys); MatFree(&Kfs); 24 } 25 else{ 26 SystemMatricesx(&Kgg, NULL, NULL, &pg,NULL, NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); 27 Reducematrixfromgtofx(&Kff,&Kfs,Kgg,femmodel->nodesets,femmodel->parameters); MatFree(&Kgg); 28 Reduceloadfromgtofx(&pf, pg, Kfs, femmodel->ys, femmodel->nodesets,femmodel->parameters);VecFree(&pg); MatFree(&Kfs); 29 } 21 30 22 31 Solverx(&uf, Kff, pf, NULL, femmodel->parameters); MatFree(&Kff); VecFree(&pf); 23 24 32 Mergesolutionfromftogx(&ug, uf,femmodel->ys,femmodel->nodesets,femmodel->parameters);VecFree(&uf); 25 33 -
issm/trunk/src/c/solvers/solver_thermal_nonlinear.cpp
r5772 r6001 34 34 int verbose=0; 35 35 bool lowmem=0; 36 bool kffpartitioning; 36 37 37 38 /*Recover parameters: */ 38 39 kflag=1; pflag=1; 39 40 femmodel->parameters->FindParam(&kffpartitioning,KffEnum); 40 41 femmodel->parameters->FindParam(&verbose,VerboseEnum); 41 42 femmodel->parameters->FindParam(&lowmem,LowmemEnum); … … 50 51 for(;;){ 51 52 52 53 SystemMatricesx(&Kgg,NULL, NULL, &pg,NULL, &melting_offset,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); 54 55 Reducematrixfromgtofx(&Kff,&Kfs,Kgg,femmodel->nodesets,femmodel->parameters); MatFree(&Kgg); 56 Reduceloadfromgtofx(&pf, pg, Kfs, femmodel->ys, femmodel->nodesets,femmodel->parameters); 57 58 VecFree(&pg); 59 MatFree(&Kfs); 53 if(kffpartitioning){ 54 SystemMatricesx(NULL,&Kff, &Kfs, NULL,&pf,&melting_offset,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); 55 Reduceloadx(pf, Kfs, femmodel->ys); MatFree(&Kfs); 56 } 57 else{ 58 SystemMatricesx(&Kgg,NULL, NULL, &pg,NULL, &melting_offset,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters); 59 Reducematrixfromgtofx(&Kff,&Kfs,Kgg,femmodel->nodesets,femmodel->parameters); MatFree(&Kgg); 60 Reduceloadfromgtofx(&pf, pg, Kfs, femmodel->ys, femmodel->nodesets,femmodel->parameters);VecFree(&pg); MatFree(&Kfs); 61 } 60 62 61 63 VecFree(&tf);
Note:
See TracChangeset
for help on using the changeset viewer.