Ignore:
Timestamp:
09/13/10 13:18:07 (15 years ago)
Author:
Eric.Larour
Message:

New solution strategy, based on creation of Kff, Kfs directly,
instead of Kgg. Idem for pf instead of pg.
Kept backwards compatibility for now, and only applied to macayeal diagnostic
in 2d. Will switch to this new strategy in the near future.

File:
1 edited

Legend:

Unmodified
Added
Removed

  • issm/trunk/src/c/solvers/solver_diagnostic_nonlinear.cpp

    r5695 r5772  
    1212void solver_diagnostic_nonlinear(FemModel* femmodel,bool conserve_loads){
    1313
    14         /*intermediary: */
    15         Mat Kgg = NULL, Kff = NULL, Kfs   = NULL;
    16         Vec ug  = NULL, uf  = NULL, old_ug= NULL, old_uf = NULL;
    17         Vec pg  = NULL, pf  = NULL;
    18         Loads* loads=NULL;
    19         int converged;
    20         int constraints_converged;
    21         int num_unstable_constraints;
    22         int count;
    23         int min_mechanical_constraints;
    24         int max_nonlinear_iterations;
    25 
    2614        /*parameters:*/
    27         int verbose=0;
     15        bool kff=false;
    2816
    2917        /*Recover parameters: */
    30         femmodel->parameters->FindParam(&verbose,VerboseEnum);
    31         femmodel->parameters->FindParam(&min_mechanical_constraints,MinMechanicalConstraintsEnum);
    32         femmodel->parameters->FindParam(&max_nonlinear_iterations,MaxNonlinearIterationsEnum);
    33        
    34         /*Were loads requested as output? : */
    35         if(conserve_loads){
    36                 loads=(Loads*)femmodel->loads->Copy(); //protect loads from being modified by the solution
    37         }
    38         else{
    39                 loads=(Loads*)femmodel->loads; //modify loads  in this solution
    40         }
     18        femmodel->parameters->FindParam(&kff,KffEnum);
    4119
    42         count=1;
    43         converged=0;
     20        if(!kff) solver_diagnostic_nonlinear_kgg(femmodel,conserve_loads);
     21        else     solver_diagnostic_nonlinear_kff(femmodel,conserve_loads);
    4422
    45         /*Start non-linear iteration using input velocity: */
    46         GetSolutionFromInputsx(&ug, femmodel->elements, femmodel->nodes, femmodel->vertices, loads, femmodel->materials, femmodel->parameters);
    47         Reducevectorgtofx(&uf, ug, femmodel->nodesets);
    48 
    49         //Update once again the solution to make sure that vx and vxold are similar (for next step in transient or steadystate)
    50         InputUpdateFromSolutionx(femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,ug);
    51 
    52         for(;;){
    53 
    54                 //save pointer to old velocity
    55                 VecFree(&old_ug);old_ug=ug;
    56                 VecFree(&old_uf);old_uf=uf;
    57 
    58                 SystemMatricesx(&Kgg, &pg,NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,loads,femmodel->materials,femmodel->parameters);
    59                
    60                 Reducematrixfromgtofx(&Kff,&Kfs,Kgg,femmodel->nodesets,femmodel->parameters); MatFree(&Kgg);
    61        
    62                 Reduceloadfromgtofx(&pf, pg, Kfs, femmodel->ys, femmodel->nodesets,femmodel->parameters); VecFree(&pg); MatFree(&Kfs);
    63 
    64                 Solverx(&uf, Kff, pf, old_uf, femmodel->parameters);
    65 
    66                 Mergesolutionfromftogx(&ug, uf,femmodel->ys,femmodel->nodesets,femmodel->parameters);
    67 
    68                 InputUpdateFromSolutionx( femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,ug);
    69 
    70                 PenaltyConstraintsx(&constraints_converged, &num_unstable_constraints, femmodel->elements,femmodel->nodes,femmodel->vertices,loads,femmodel->materials,femmodel->parameters);
    71                 if(verbose)_printf_("   number of unstable constraints: %i\n",num_unstable_constraints);
    72 
    73                 convergence(&converged,Kff,pf,uf,old_uf,femmodel->parameters); MatFree(&Kff);VecFree(&pf);
    74                
    75                 InputUpdateFromConstantx( femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,converged,ConvergedEnum);
    76 
    77                 //rift convergence
    78                 if (!constraints_converged) {
    79                         if (converged){
    80                                 if (num_unstable_constraints <= min_mechanical_constraints) converged=1;
    81                                 else converged=0;
    82                         }
    83                 }
    84 
    85                 /*Increase count: */
    86                 count++;
    87                 if(converged==1)break;
    88                 if(count>=max_nonlinear_iterations){
    89                         _printf_("   maximum number of iterations (%i) exceeded\n",max_nonlinear_iterations);
    90                         break;
    91                 }
    92         }
    93 
    94         /*clean-up*/
    95         if(conserve_loads) delete loads;
    96         VecFree(&uf);
    97         VecFree(&ug);
    98         VecFree(&old_uf);
    99         VecFree(&old_ug);
    100        
    10123}
Note: See TracChangeset for help on using the changeset viewer.