Changeset 11679


Ignore:
Timestamp:
03/12/12 09:35:33 (13 years ago)
Author:
Eric.Larour
Message:

Large commit, to strip out Petsc from the code. Created new Matrix and Vector
wrapper objects to the Petsc Mat and Vec objects. These wrappers make it possible
to hook up a different package to handle matrix and vector operations.

Location:
issm/trunk-jpl/src
Files:
1 added
73 edited

Legend:

Unmodified
Added
Removed

  • issm/trunk-jpl/src/c/Makefile.am

    r11670 r11679  
    344344                                        ./modules/ResetCoordinateSystemx/ResetCoordinateSystemx.cpp\
    345345                                        ./modules/Solverx/Solverx.cpp\
     346                                        ./modules/Solverx/SolverxPetsc.cpp\
    346347                                        ./modules/Solverx/Solverx.h\
    347348                                        ./modules/Solverx/DofTypesToIndexSet.cpp\

  • issm/trunk-jpl/src/c/modules/CreateJacobianMatrixx/CreateJacobianMatrixx.cpp

    r11327 r11679  
    1010#include "../../EnumDefinitions/EnumDefinitions.h"
    1111
    12 void CreateJacobianMatrixx(Mat* pJff,Elements* elements,Nodes* nodes, Vertices* vertices,Loads* loads,Materials* materials, Parameters* parameters,double kmax){
     12void CreateJacobianMatrixx(Matrix** pJff,Elements* elements,Nodes* nodes, Vertices* vertices,Loads* loads,Materials* materials, Parameters* parameters,double kmax){
    1313       
    1414        int      i,connectivity;
    … …  
    1717        Element *element = NULL;
    1818        Load    *load    = NULL;
    19         Mat      Jff     = NULL;
     19        Matrix*  Jff     = NULL;
    2020
    2121        /*Checks*/
    … …  
    2929
    3030        /*Initialize Jacobian Matrix*/
    31         Jff=NewMat(fsize,fsize,connectivity,numberofdofspernode);
     31        Jff=new Matrix(fsize,fsize,connectivity,numberofdofspernode);
    3232       
    3333        /*Create and assemble matrix*/
    … …  
    4141                if(load->InAnalysis(configuration_type)) load->PenaltyCreateJacobianMatrix(Jff,kmax);
    4242        }
    43         MatAssemblyBegin(Jff,MAT_FINAL_ASSEMBLY);
    44         MatAssemblyEnd(Jff,MAT_FINAL_ASSEMBLY);
     43        Jff->Assemble();
    4544
    4645        /*Assign output pointer*/

  • issm/trunk-jpl/src/c/modules/CreateJacobianMatrixx/CreateJacobianMatrixx.h

    r11327 r11679  
    1010
    1111/* local prototypes: */
    12 void CreateJacobianMatrixx(Mat* pJff,Elements* elements,Nodes* nodes, Vertices* vertices,Loads* loads,Materials* materials, Parameters* parameters,double kmax);
     12void CreateJacobianMatrixx(Matrix** pJff,Elements* elements,Nodes* nodes, Vertices* vertices,Loads* loads,Materials* materials, Parameters* parameters,double kmax);
    1313
    1414#endif  /* _CREATEJACOBIANMATRIXX_H */

  • issm/trunk-jpl/src/c/modules/CreateNodalConstraintsx/CreateNodalConstraintsx.cpp

    r8816 r11679  
    1010#include "../../EnumDefinitions/EnumDefinitions.h"
    1111
    12 void CreateNodalConstraintsx( Vec* pys, Nodes* nodes,int configuration_type){
     12void CreateNodalConstraintsx( Vector** pys, Nodes* nodes,int configuration_type){
    1313
    1414        int i;
    … …  
    1818
    1919        /*output: */
    20         Vec ys=NULL;
     20        Vector* ys=NULL;
    2121
    2222        /*figure out how many dofs we have: */
    … …  
    2424
    2525        /*allocate:*/
    26         ys=NewVec(numberofdofs);
     26        ys=new Vector(numberofdofs);
    2727
    2828        /*go through all nodes, and for the ones corresponding to this configuration_type, fill the
    … …  
    3636
    3737        /*Assemble: */
    38         VecAssemblyBegin(ys);
    39         VecAssemblyEnd(ys);
     38        ys->Assemble();
    4039
    4140        /*Assign output pointers: */

  • issm/trunk-jpl/src/c/modules/CreateNodalConstraintsx/CreateNodalConstraintsx.h

    r8816 r11679  
    99
    1010/* local prototypes: */
    11 void CreateNodalConstraintsx( Vec* pys, Nodes* nodes,int configuration_type);
     11void CreateNodalConstraintsx( Vector** pys, Nodes* nodes,int configuration_type);
    1212
    1313#endif  /* _CREATENODALCONSTRAINTSX_H */

  • issm/trunk-jpl/src/c/modules/GetSolutionFromInputsx/GetSolutionFromInputsx.cpp

    r8224 r11679  
    99#include "../../EnumDefinitions/EnumDefinitions.h"
    1010
    11 void    GetSolutionFromInputsx( Vec* psolution, Elements* elements,Nodes* nodes, Vertices* vertices,Loads* loads, Materials* materials,  Parameters* parameters){
     11void    GetSolutionFromInputsx( Vector** psolution, Elements* elements,Nodes* nodes, Vertices* vertices,Loads* loads, Materials* materials,  Parameters* parameters){
    1212
    1313        /*intermediary: */
    … …  
    1919
    2020        /*output: */
    21         Vec solution=NULL;
     21        Vector* solution=NULL;
    2222
    2323        /*retrive parameters: */
    … …  
    2929       
    3030        /*Initialize solution: */
    31         solution=NewVec(gsize);
     31        solution=new Vector(gsize);
    3232       
    3333        /*Go through elements and plug solution: */
    … …  
    3838
    3939        /*Assemble vector: */
    40         VecAssemblyBegin(solution);
    41         VecAssemblyEnd(solution);
     40        solution->Assemble();
    4241
    4342        /*Assign output pointers:*/

  • issm/trunk-jpl/src/c/modules/GetSolutionFromInputsx/GetSolutionFromInputsx.h

    r4236 r11679  
    1010
    1111/* local prototypes: */
    12 void GetSolutionFromInputsx( Vec* psolution, Elements* elements,Nodes* nodes, Vertices* vertices,Loads* loads, Materials* materials,  Parameters* parameters);
     12void GetSolutionFromInputsx( Vector** psolution, Elements* elements,Nodes* nodes, Vertices* vertices,Loads* loads, Materials* materials,  Parameters* parameters);
    1313
    1414#endif  /* _GETSOLUTIONFROMINPUTSXX_H */

  • issm/trunk-jpl/src/c/modules/InputUpdateFromSolutionx/InputUpdateFromSolutionx.cpp

    r4573 r11679  
    99#include "../../EnumDefinitions/EnumDefinitions.h"
    1010
    11 void InputUpdateFromSolutionx( Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials, Parameters* parameters,Vec solution){
     11void InputUpdateFromSolutionx( Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials, Parameters* parameters,Vector* solution){
    1212
    1313        double* serial_solution=NULL;
    1414
    1515        /*Serialize solution, so that elements can index into it on every CPU: */
    16         VecToMPISerial(&serial_solution,solution);
     16        serial_solution=solution->ToMPISerial();
    1717
    1818        /*Call overloaded form of InputUpdateFromSolutionx: */

  • issm/trunk-jpl/src/c/modules/InputUpdateFromSolutionx/InputUpdateFromSolutionx.h

    r4236 r11679  
    1010
    1111/* local prototypes: */
    12 void            InputUpdateFromSolutionx( Elements* elements,Nodes* nodes, Vertices* vertices,Loads* loads, Materials* materials,  Parameters* parameters,Vec solution);
     12void            InputUpdateFromSolutionx( Elements* elements,Nodes* nodes, Vertices* vertices,Loads* loads, Materials* materials,  Parameters* parameters,Vector* solution);
    1313void        InputUpdateFromSolutionx( Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials, Parameters* parameters,double* solution);
    1414
    1515//with timestep
    16 void            InputUpdateFromSolutionx( Elements* elements,Nodes* nodes, Vertices* vertices,Loads* loads, Materials* materials,  Parameters* parameters,Vec solution,int timestep);
     16void            InputUpdateFromSolutionx( Elements* elements,Nodes* nodes, Vertices* vertices,Loads* loads, Materials* materials,  Parameters* parameters,Vector* solution,int timestep);
    1717void        InputUpdateFromSolutionx( Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials, Parameters* parameters,double* solution, int timestep);
    1818

  • issm/trunk-jpl/src/c/modules/Mergesolutionfromftogx/Mergesolutionfromftogx.cpp

    r9761 r11679  
    77#include "./Mergesolutionfromftogx.h"
    88
    9 void    Mergesolutionfromftogx( Vec* pug, Vec uf, Vec ys, Nodes* nodes, Parameters* parameters, bool flag_ys0){
     9void    Mergesolutionfromftogx( Vector** pug, Vector* uf, Vector* ys, Nodes* nodes, Parameters* parameters, bool flag_ys0){
    1010
    1111        /*output: */
    12         Vec ug=NULL;
     12        Vector* ug=NULL;
    1313
    1414        /*intermediary: */
    … …  
    2929        if(ssize){
    3030                if(flag_ys0){
    31                         VecSet(ys,0.0);
     31                        ys->Set(0.0);
    3232                }
    3333        }
    3434
    3535        /*initialize ug: */
    36         ug=NewVec(gsize);
     36        ug=new Vector(gsize);
    3737
    3838        /*Merge f set back into g set: */

  • issm/trunk-jpl/src/c/modules/Mergesolutionfromftogx/Mergesolutionfromftogx.h

    r8799 r11679  
    99
    1010/* local prototypes: */
    11 void    Mergesolutionfromftogx( Vec* pug, Vec uf, Vec ys, Nodes* nodes, Parameters* parameters, bool flag_ys0=false);
     11void    Mergesolutionfromftogx( Vector** pug, Vector* uf, Vector* ys, Nodes* nodes, Parameters* parameters, bool flag_ys0=false);
    1212
    1313#endif  /* _MERGESOLUTIONFROMFTOGX_H */

  • issm/trunk-jpl/src/c/modules/Reduceloadx/Reduceloadx.cpp

    r9761 r11679  
    1212#include "../../io/io.h"
    1313
    14 void    Reduceloadx( Vec pf, Mat Kfs, Vec y_s,bool flag_ys0){
     14void    Reduceloadx( Vector* pf, Matrix* Kfs, Vector* y_s,bool flag_ys0){
    1515
    1616        /*intermediary*/
    17         Vec         y_s0   = NULL;
    18         Vec         Kfsy_s = NULL;
     17        Vector*     y_s0   = NULL;
     18        Vector*     Kfsy_s = NULL;
    1919        int         Kfsm,Kfsn;
    2020        int         global_m,global_n;
    21         PetscScalar a;
    2221        bool        fromlocalsize = true;
    2322        int         verbose;
    … …  
    2524        _printf_(VerboseModule(),"   Dirichlet lifting applied to load vector\n");
    2625
    27         MatGetSize(Kfs,&global_m,&global_n);
     26        Kfs->GetSize(&global_m,&global_n);
    2827        if(pf && global_m*global_n){
    2928
    … …  
    3231
    3332                /*pf = pf - Kfs * y_s;*/
    34                 MatGetLocalSize(Kfs,&Kfsm,&Kfsn);
    35                 Kfsy_s=NewVec(Kfsm,fromlocalsize);
     33                Kfs->GetLocalSize(&Kfsm,&Kfsn);
     34                Kfsy_s=new Vector(Kfsm,fromlocalsize);
    3635                if (flag_ys0){
    3736
    3837                        /*Create y_s0, full of 0: */
    39                         VecDuplicate(y_s,&y_s0);
    40                         VecSet(y_s0,0.0);
    41                         VecAssemblyBegin(y_s0);
    42                         VecAssemblyEnd(y_s0);
     38                        y_s0=y_s->Duplicate();
     39                        y_s0->Set(0.0);
     40                        y_s0->Assemble();
    4341
    44                         MatMultPatch(Kfs,y_s0,Kfsy_s);
     42                        Kfs->MatMult(y_s0,Kfsy_s);
    4543                }
    4644                else{
    47                         MatMultPatch(Kfs,y_s,Kfsy_s);
     45                        Kfs->MatMult(y_s,Kfsy_s);
    4846                }
    4947
    50                 a=-1;
    51                 VecAXPY(pf,a,Kfsy_s); 
     48                pf->AXPY(Kfsy_s,-1);
    5249        }
    5350
    5451
    5552        /*Free ressources and return*/
    56         VecFree(&y_s0);
    57         VecFree(&Kfsy_s);
     53        delete y_s0;
     54        delete Kfsy_s;
    5855
    5956}

  • issm/trunk-jpl/src/c/modules/Reduceloadx/Reduceloadx.h

    r6580 r11679  
    99
    1010/* local prototypes: */
    11 void    Reduceloadx( Vec pf, Mat Kfs, Vec ys,bool flag_ys0=false);
     11void    Reduceloadx( Vector* pf, Matrix* Kfs, Vector* ys,bool flag_ys0=false);
    1212
    1313#endif  /* _REDUCELOADX_H */

  • issm/trunk-jpl/src/c/modules/Reducevectorgtofx/Reducevectorgtofx.cpp

    r8817 r11679  
    66#include "./Reducevectorgtofx.h"
    77 
    8 void Reducevectorgtofx(Vec* puf, Vec ug, Nodes* nodes,Parameters* parameters){
     8void Reducevectorgtofx(Vector** puf, Vector* ug, Nodes* nodes,Parameters* parameters){
    99
    1010        /*output: */
    11         Vec uf=NULL;
     11        Vector* uf=NULL;
    1212
    1313        /*variables: */
    … …  
    2626        else{
    2727                /*allocate: */
    28                 uf=NewVec(fsize);
     28                uf=new Vector(fsize);
    2929
    3030                if(nodes->NumberOfNodes(configuration_type)){
    3131
    3232                        /*serialize ug, so nodes can index into it: */
    33                         VecToMPISerial(&ug_serial,ug);
     33                        ug_serial=ug->ToMPISerial();
    3434
    3535                        /*Go through all nodes, and ask them to retrieve values from ug, and plug them into uf: */
    … …  
    4747                }
    4848                /*Assemble vector: */
    49                 VecAssemblyBegin(uf);
    50                 VecAssemblyEnd(uf);
     49                uf->Assemble();
    5150        }
    5251

  • issm/trunk-jpl/src/c/modules/Reducevectorgtofx/Reducevectorgtofx.h

    r8799 r11679  
    1010
    1111/* local prototypes: */
    12 void Reducevectorgtofx(Vec* puf, Vec ug, Nodes* nodes,Parameters* parameters);
     12void Reducevectorgtofx(Vector** puf, Vector* ug, Nodes* nodes,Parameters* parameters);
    1313
    1414#endif  /* _REDUCEVECTORGTOFX_H */

  • issm/trunk-jpl/src/c/modules/Reducevectorgtosx/Reducevectorgtosx.cpp

    r8816 r11679  
    66#include "./Reducevectorgtosx.h"
    77
    8 void Reducevectorgtosx(Vec* pys, Vec yg, Nodes* nodes,Parameters* parameters){
     8void Reducevectorgtosx(Vector** pys, Vector* yg, Nodes* nodes,Parameters* parameters){
    99
    1010        /*output: */
    11         Vec ys=NULL;
     11        Vector* ys=NULL;
    1212
    1313        /*variables: */
    … …  
    2626        else{
    2727                /*allocate: */
    28                 ys=NewVec(ssize);
     28                ys=new Vector(ssize);
    2929
    3030                if(nodes->NumberOfNodes(configuration_type)){
    3131
    3232                        /*serialize yg, so nodes can index into it: */
    33                         VecToMPISerial(&yg_serial,yg);
     33                        yg_serial=yg->ToMPISerial();
    3434
    3535                        /*Go throygh all nodes, and ask them to retrieve values from yg, and plyg them into ys: */
    … …  
    4747                }
    4848                /*Assemble vector: */
    49                 VecAssemblyBegin(ys);
    50                 VecAssemblyEnd(ys);
     49                ys->Assemble();
    5150        }
    5251

  • issm/trunk-jpl/src/c/modules/Reducevectorgtosx/Reducevectorgtosx.h

    r8799 r11679  
    1010
    1111/* local prototypes: */
    12 void Reducevectorgtosx(Vec* pys, Vec yg, Nodes* nodes,Parameters* parameters);
     12void Reducevectorgtosx(Vector** pys, Vector* yg, Nodes* nodes,Parameters* parameters);
    1313
    1414#endif  /* _REDUCEVECTORGTOSX_H */

  • issm/trunk-jpl/src/c/modules/Solverx/Solverx.cpp

    r11469 r11679  
    1414#endif
    1515
    16 void    Solverx(Vec* puf, Mat Kff, Vec pf, Vec uf0,Vec df, Parameters* parameters){
     16void    Solverx(Vector** puf, Matrix* Kff, Vector* pf, Vector* uf0,Vector* df, Parameters* parameters){
    1717
    18         /*Output: */
    19         Vec        uf               = NULL;
     18        /*output: */
     19        Vector* uf=NULL;
     20        uf=new Vector();
    2021
    21         /*Intermediary: */
    22         int        local_m,local_n,global_m,global_n;
    23         int        analysis_type;
     22        #ifdef _HAVE_PETSC_
     23        Vec uf0_vector=NULL;
     24        if (uf0)uf0_vector=uf0->vector;
    2425
    25         /*Solver */
    26         KSP        ksp              = NULL;
    27         PC         pc               = NULL;
    28         int        iteration_number;
    29         int        solver_type;
    30         bool       fromlocalsize    = true;
    31         #if _PETSC_MAJOR_ < 3 || (_PETSC_MAJOR_ == 3 && _PETSC_MINOR_ < 2)
    32         PetscTruth flag,flg;
     26        SolverxPetsc(&uf->vector, Kff->matrix, pf->vector, uf0_vector, df->vector, parameters);
     27        VecGetSize(uf->vector,&uf->M);
    3328        #else
    34         PetscBool flag,flg;
     29        _error_("not supported yet!");
    3530        #endif
    3631
    37         /*Stokes: */
    38         IS         isv=NULL;
    39         IS         isp=NULL;
    40 
    41         #if _PETSC_MAJOR_ >= 3
    42         char ksp_type[50];
    43         #endif
    44 
    45 
    46         /*Display message*/
    47         _printf_(VerboseModule(),"   Solving\n");
    48         #if _PETSC_MAJOR_ < 3 || (_PETSC_MAJOR_ == 3 && _PETSC_MINOR_ < 2)
    49         if(VerboseSolver())PetscOptionsPrint(stdout);
    50         #else
    51         if(VerboseSolver())PetscOptionsView(PETSC_VIEWER_STDOUT_WORLD);
    52         #endif
    53 
    54         /*First, check that f-set is not NULL, i.e. model is fully constrained: {{{*/
    55         _assert_(Kff);
    56         MatGetSize(Kff,&global_m,&global_n); _assert_(global_m==global_m);
    57         if(!global_n){
    58                 *puf=NULL; return;
    59         }
    60         /*}}}*/
    61         /*Initial guess logic here: {{{1*/
    62         /*Now, check that we are not giving an initial guess to the solver, if we are running a direct solver: */
    63         #if _PETSC_MAJOR_ >= 3
    64         PetscOptionsGetString(PETSC_NULL,"-ksp_type",ksp_type,49,&flg);
    65         if (strcmp(ksp_type,"preonly")==0)uf0=NULL;
    66         #endif
    67 
    68         /*If initial guess for the solution exists, use it to create uf, otherwise,
    69          * duplicate the right hand side so that the solution vector has the same structure*/
    70         if(uf0){
    71                 VecDuplicate(uf0,&uf); VecCopy(uf0,uf);
    72         }
    73         else{
    74                 MatGetLocalSize(Kff,&local_m,&local_n);uf=NewVec(local_n,fromlocalsize);
    75         }
    76         /*}}}*/
    77         /*Process petsc options to see if we are using special types of external solvers: {{{1*/
    78         PetscOptionsDetermineSolverType(&solver_type);
    79 
    80         /*In serial mode, the matrices have been loaded as MPIAIJ or AIJ matrices.
    81          We need to convert them if we are going to run the solvers successfully: */
    82         #ifdef _SERIAL_
    83         #if _PETSC_MAJOR_ == 2
    84         if (solver_type==MUMPSPACKAGE_LU){
    85                 /*Convert Kff to MATTAIJMUMPS: */
    86                 MatConvert(Kff,MATAIJMUMPS,MAT_REUSE_MATRIX,&Kff);
    87         }
    88         if (solver_type==MUMPSPACKAGE_CHOL){
    89                 /*Convert Kff to MATTSBAIJMUMPS: */
    90                 MatConvert(Kff,MATSBAIJMUMPS,MAT_REUSE_MATRIX,&Kff);
    91         }
    92         if (solver_type==SPOOLESPACKAGE_LU){
    93                 /*Convert Kff to MATTSBAIJMUMPS: */
    94                 MatConvert(Kff,MATAIJSPOOLES,MAT_REUSE_MATRIX,&Kff);
    95         }
    96         if (solver_type==SPOOLESPACKAGE_CHOL){
    97                 /*Convert Kff to MATTSBAIJMUMPS: */
    98                 MatConvert(Kff,MATSBAIJSPOOLES,MAT_REUSE_MATRIX,&Kff);
    99         }
    100         if (solver_type==SUPERLUDISTPACKAGE){
    101                 /*Convert Kff to MATTSBAIJMUMPS: */
    102                 MatConvert(Kff,MATSUPERLU_DIST,MAT_REUSE_MATRIX,&Kff);
    103         }
    104         if (solver_type==StokesSolverEnum){
    105                 _error_("Petsc 2 does not support multi-physics solvers");
    106         }
    107         #endif
    108         #endif
    109         /*}}}*/
    110         /*Check the solver is available: {{{1*/
    111         if(solver_type==MUMPSPACKAGE_LU || solver_type==MUMPSPACKAGE_CHOL){
    112         #if _PETSC_MAJOR_ >=3
    113                 #ifndef _HAVE_MUMPS_
    114                 _error_("requested MUMPS solver, which was not compiled into ISSM!\n");
    115                 #endif
    116 
    117         #endif
    118         }
    119         /*}}}*/
    120         /*Prepare solver:{{{1*/
    121         KSPCreate(MPI_COMM_WORLD,&ksp);
    122         KSPSetOperators(ksp,Kff,Kff,DIFFERENT_NONZERO_PATTERN);
    123         KSPSetFromOptions(ksp);
    124 
    125         #if defined(_SERIAL_) && _PETSC_MAJOR_==3
    126         /*Specific solver?: */
    127         KSPGetPC(ksp,&pc);
    128         if (solver_type==MUMPSPACKAGE_LU){
    129                 #if _PETSC_MINOR_==1
    130                 PCFactorSetMatSolverPackage(pc,MAT_SOLVER_MUMPS);
    131                 #else
    132                 PCFactorSetMatSolverPackage(pc,MATSOLVERMUMPS);
    133                 #endif
    134         }
    135         #endif
    136 
    137         #if defined(_PARALLEL_) && _PETSC_MAJOR_==3
    138         /*Stokes: */
    139         if (solver_type==StokesSolverEnum){
    140                 /*Make indices out of doftypes: */
    141                 if(!df)_error_("need doftypes for Stokes solver!\n");
    142                 DofTypesToIndexSet(&isv,&isp,df,StokesSolverEnum);
    143 
    144                 /*Set field splits: */
    145                 KSPGetPC(ksp,&pc);
    146                 #if _PETSC_MINOR_==1
    147                 PCFieldSplitSetIS(pc,isv);
    148                 PCFieldSplitSetIS(pc,isp);
    149                 #else
    150                 PCFieldSplitSetIS(pc,PETSC_NULL,isv);
    151                 PCFieldSplitSetIS(pc,PETSC_NULL,isp);
    152                 #endif
    153 
    154         }
    155         #endif
    156 
    157         /*}}}*/
    158         /*If there is an initial guess for the solution, use it, except if we are using the MUMPS direct solver, where any initial guess will crash Petsc: {{{1*/
    159         if (uf0){
    160                 if( (solver_type!=MUMPSPACKAGE_LU) && (solver_type!=MUMPSPACKAGE_CHOL) && (solver_type!=SPOOLESPACKAGE_LU) && (solver_type!=SPOOLESPACKAGE_CHOL) && (solver_type!=SUPERLUDISTPACKAGE)){
    161                         KSPSetInitialGuessNonzero(ksp,PETSC_TRUE);
    162                 }
    163         }
    164         /*}}}*/
    165        
    166         if(VerboseSolver())KSPView(ksp,PETSC_VIEWER_STDOUT_WORLD);
    167 
    168         /*Solve: */
    169         KSPSolve(ksp,pf,uf);
    170        
    171         /*Check convergence*/
    172         KSPGetIterationNumber(ksp,&iteration_number);
    173         if (iteration_number<0) _error_("%s%i"," Solver diverged at iteration number: ",-iteration_number);
    174 
    175         /*Free resources:*/
    176         KSPFree(&ksp);
    177                
    178         /*Assign output pointers:*/
     32        /*Assign output pointers: */
    17933        *puf=uf;
    18034}

  • issm/trunk-jpl/src/c/modules/Solverx/Solverx.h

    r7391 r11679  
    99
    1010/* local prototypes: */
    11 void    Solverx( Vec* puf, Mat Kff, Vec pf, Vec uf0,Vec df,Parameters* parameters);
     11void    Solverx(Vector** puf, Matrix* Kff, Vector* pf, Vector* uf0,Vector* df, Parameters* parameters);
     12void    SolverxPetsc(Vec* puf, Mat Kff, Vec pf, Vec uf0,Vec df, Parameters* parameters);
    1213void    DofTypesToIndexSet(IS* pisv, IS* pisp, Vec df,int typeenum);
    1314#endif  /* _SOLVERX_H */

  • issm/trunk-jpl/src/c/modules/SystemMatricesx/SystemMatricesx.cpp

    r11324 r11679  
    1010#include "../../EnumDefinitions/EnumDefinitions.h"
    1111
    12 void SystemMatricesx(Mat* pKff, Mat* pKfs, Vec* ppf, Vec* pdf, double* pkmax,Elements* elements,Nodes* nodes, Vertices* vertices,Loads* loads,Materials* materials, Parameters* parameters,bool kflag,bool pflag,bool penalty_kflag,bool penalty_pflag){
     12void SystemMatricesx(Matrix** pKff, Matrix** pKfs, Vector** ppf, Vector** pdf, double* pkmax,Elements* elements,Nodes* nodes, Vertices* vertices,Loads* loads,Materials* materials, Parameters* parameters,bool kflag,bool pflag,bool penalty_kflag,bool penalty_pflag){
    1313       
    1414        /*intermediary: */
    … …  
    2121       
    2222        /*output: */
    23         Mat    Kff  = NULL;
    24         Mat    Kfs  = NULL;
    25         Vec    pf   = NULL;
    26         Vec    df=NULL;
     23        Matrix*    Kff  = NULL;
     24        Matrix*    Kfs  = NULL;
     25        Vector*    pf   = NULL;
     26        Vector*    df=NULL;
    2727        double kmax = 0;
    2828
    … …  
    4949        if(kflag){
    5050
    51                 Kff=NewMat(fsize,fsize,connectivity,numberofdofspernode);
    52                 Kfs=NewMat(fsize,ssize,connectivity,numberofdofspernode);
    53                 df=NewVec(fsize);
     51                Kff=new Matrix(fsize,fsize,connectivity,numberofdofspernode);
     52                Kfs=new Matrix(fsize,ssize,connectivity,numberofdofspernode);
     53                df=new Vector(fsize);
    5454
    5555                /*Fill stiffness matrix from elements: */
    … …  
    6666
    6767                /*Assemble matrix and doftypes and compress matrix to save memory: */
    68                 MatAssemblyBegin(Kff,MAT_FINAL_ASSEMBLY);
    69                 MatAssemblyEnd(Kff,MAT_FINAL_ASSEMBLY);
    70                 #if _PETSC_MAJOR_ == 2
    71                 MatCompress(Kff);
    72                 #endif
     68                Kff->Assemble();
     69                Kfs->Assemble();
     70                df->Assemble();
    7371
    74                 MatAssemblyBegin(Kfs,MAT_FINAL_ASSEMBLY);
    75                 MatAssemblyEnd(Kfs,MAT_FINAL_ASSEMBLY);
    76                 #if _PETSC_MAJOR_ == 2
    77                 MatCompress(Kfs);
    78                 #endif
    79                 VecAssemblyBegin(df);
    80                 VecAssemblyEnd(df);
    8172               
    8273        }
    … …  
    8475        if(pflag){
    8576
    86                 pf=NewVec(fsize);
     77                pf=new Vector(fsize);
    8778
    8879                /*Fill right hand side vector, from elements: */
    … …  
    9889                }
    9990
    100                 VecAssemblyBegin(pf);
    101                 VecAssemblyEnd(pf);
     91                pf->Assemble();
    10292        }
    10393
    10494        /*Now, figure out maximum value of K_gg, so that we can penalize it correctly: */
    105         MatNorm(Kff,NORM_INFINITY,&kmax);
     95        kmax=Kff->Norm(NORM_INFINITY);
    10696
    10797        /*Now, deal with penalties*/
    … …  
    115105
    116106                /*Assemble matrix and compress matrix to save memory: */
    117                 MatAssemblyBegin(Kff,MAT_FINAL_ASSEMBLY);
    118                 MatAssemblyEnd(Kff,MAT_FINAL_ASSEMBLY);
    119                 #if _PETSC_MAJOR_ == 2
    120                 MatCompress(Kff);
    121                 #endif
    122 
    123                 MatAssemblyBegin(Kfs,MAT_FINAL_ASSEMBLY);
    124                 MatAssemblyEnd(Kfs,MAT_FINAL_ASSEMBLY);
    125                 #if _PETSC_MAJOR_ == 2
    126                 MatCompress(Kfs);
    127                 #endif
     107                Kff->Assemble();
     108                Kfs->Assemble();
    128109        }
    129110
    … …  
    137118                }
    138119
    139                 VecAssemblyBegin(pf);
    140                 VecAssemblyEnd(pf);
     120                pf->Assemble();
     121                pf->Assemble();
    141122        }
    142123
    143124        /*Assign output pointers: */
    144125        if(pKff) *pKff=Kff;
    145         else      MatFree(&Kff);
     126        else      delete Kff;
    146127        if(pKfs) *pKfs=Kfs;
    147         else      MatFree(&Kfs);
     128        else      delete Kfs;
    148129        if(ppf)  *ppf=pf;
    149         else      VecFree(&pf);
     130        else      delete pf;
    150131        if(pdf)  *pdf=df;
    151         else      VecFree(&df);
     132        else      delete df;
    152133        if(pkmax) *pkmax=kmax;
    153134}

  • issm/trunk-jpl/src/c/modules/SystemMatricesx/SystemMatricesx.h

    r8799 r11679  
    1010
    1111/* local prototypes: */
    12 void SystemMatricesx(Mat* pKff, Mat* pKfs, Vec* ppf, Vec* pdf, double* pkmax,Elements* elements,Nodes* nodes, Vertices* vertices,Loads* loads,Materials* materials, Parameters* parameters,
     12void SystemMatricesx(Matrix** pKff, Matrix** pKfs, Vector** ppf, Vector** pdf, double* pkmax,Elements* elements,Nodes* nodes, Vertices* vertices,Loads* loads,Materials* materials, Parameters* parameters,
    1313                        bool kflag=true,bool pflag=true,bool penalty_kflag=true,bool penalty_pflag=true);
    1414

  • issm/trunk-jpl/src/c/modules/UpdateDynamicConstraintsx/UpdateDynamicConstraintsx.cpp

    r9883 r11679  
    1010#include "../../EnumDefinitions/EnumDefinitions.h"
    1111
    12 void UpdateDynamicConstraintsx(Constraints* constraints,Nodes* nodes,Parameters* parameters,Vec yg){
     12void UpdateDynamicConstraintsx(Constraints* constraints,Nodes* nodes,Parameters* parameters,Vector* yg){
    1313       
    1414        int configuration_type;
    … …  
    1919
    2020        /*serialize yg, so nodes can index into it: */
    21         VecToMPISerial(&yg_serial,yg);
     21        yg_serial=yg->ToMPISerial();
    2222
    2323        for(int i=0;i<constraints->Size();i++){

  • issm/trunk-jpl/src/c/modules/UpdateDynamicConstraintsx/UpdateDynamicConstraintsx.h

    r9298 r11679  
    99#include "../../objects/objects.h"
    1010
    11 void UpdateDynamicConstraintsx(Constraints* constraints,Nodes* nodes,Parameters* parameters,Vec yg);
     11void UpdateDynamicConstraintsx(Constraints* constraints,Nodes* nodes,Parameters* parameters,Vector* yg);
    1212
    1313#endif  /* _UPDATESPCSX_H */

  • issm/trunk-jpl/src/c/modules/VecMergex/VecMergex.cpp

    r8809 r11679  
    1010#include "../../EnumDefinitions/EnumDefinitions.h"
    1111
    12 void VecMergex(Vec ug, Vec uf, Nodes* nodes, Parameters* parameters, int SetEnum){
     12void VecMergex(Vector* ug, Vector* uf, Nodes* nodes, Parameters* parameters, int SetEnum){
    1313
    1414        /*variables: */
    … …  
    2121       
    2222        /*serialize uf: */
    23         VecToMPISerial(&uf_serial,uf);
     23        uf_serial=uf->ToMPISerial();
     24
    2425
    2526        /*Do we have any nodes for this configuration? :*/
    … …  
    4344
    4445        /*Assemble vector: */
    45         VecAssemblyBegin(ug);
    46         VecAssemblyEnd(ug);
    47 
     46        ug->Assemble();
    4847}

  • issm/trunk-jpl/src/c/modules/VecMergex/VecMergex.h

    r8799 r11679  
    1010
    1111/* local prototypes: */
    12 void VecMergex(Vec ug, Vec uf, Nodes* nodes, Parameters* parameters, int SetEnum);
     12void VecMergex(Vector* ug, Vector* uf, Nodes* nodes, Parameters* parameters, int SetEnum);
    1313
    1414#endif  /* _VECMERGEX_H */

  • issm/trunk-jpl/src/c/objects/Elements/Element.h

    r11656 r11679  
    1616class Parameters;
    1717class Patch;
     18class Matrix;
     19class Vector;
    1820
    1921#include "../../toolkits/toolkits.h"
    … …  
    2830                virtual void   Configure(Elements* elements,Loads* loads,DataSet* nodes,Materials* materials,Parameters* parameters)=0;
    2931                virtual void   SetCurrentConfiguration(Elements* elements,Loads* loads,DataSet* nodes,Materials* materials,Parameters* parameters)=0;
    30                 virtual void   CreateKMatrix(Mat Kff, Mat Kfs,Vec df)=0;
    31                 virtual void   CreatePVector(Vec pf)=0;
    32                 virtual void   CreateJacobianMatrix(Mat Jff)=0;
    33                 virtual void   GetSolutionFromInputs(Vec solution)=0;
     32                virtual void   CreateKMatrix(Matrix* Kff, Matrix*  Kfs,Vector* df)=0;
     33                virtual void   CreatePVector(Vector* pf)=0;
     34                virtual void   CreateJacobianMatrix(Matrix* Jff)=0;
     35                virtual void   GetSolutionFromInputs(Vector* solution)=0;
    3436                virtual int    GetNodeIndex(Node* node)=0;
    3537                virtual int    Sid()=0;

  • issm/trunk-jpl/src/c/objects/Elements/Penta.cpp

    r11656 r11679  
    567567/*}}}*/
    568568/*FUNCTION Penta::CreateKMatrix {{{1*/
    569 void  Penta::CreateKMatrix(Mat Kff, Mat Kfs,Vec df){
     569void  Penta::CreateKMatrix(Matrix* Kff, Matrix* Kfs,Vector* df){
    570570
    571571        /*retrieve parameters: */
    … …  
    671671/*}}}*/
    672672/*FUNCTION Penta::CreatePVector {{{1*/
    673 void  Penta::CreatePVector(Vec pf){
     673void  Penta::CreatePVector(Vector* pf){
    674674
    675675        /*retrive parameters: */
    … …  
    773773/*}}}*/
    774774/*FUNCTION Penta::CreateJacobianMatrix{{{1*/
    775 void  Penta::CreateJacobianMatrix(Mat Jff){
     775void  Penta::CreateJacobianMatrix(Matrix* Jff){
    776776
    777777        /*retrieve parameters: */
    … …  
    10911091/*}}}*/
    10921092/*FUNCTION Penta::GetSolutionFromInputs{{{1*/
    1093 void  Penta::GetSolutionFromInputs(Vec solution){
     1093void  Penta::GetSolutionFromInputs(Vector* solution){
    10941094
    10951095        int analysis_type;
    … …  
    41804180/*}}}*/
    41814181/*FUNCTION Penta::GetSolutionFromInputsThermal{{{1*/
    4182 void  Penta::GetSolutionFromInputsThermal(Vec solution){
     4182void  Penta::GetSolutionFromInputsThermal(Vector* solution){
    41834183
    41844184        const int    numdof=NDOF1*NUMVERTICES;
    … …  
    42034203
    42044204        /*Add value to global vector*/
    4205         VecSetValues(solution,numdof,doflist,(const double*)values,INSERT_VALUES);
     4205        solution->SetValues(numdof,doflist,values,INSERT_VALUES);
    42064206
    42074207        /*Free ressources:*/
    … …  
    42114211/*}}}*/
    42124212/*FUNCTION Penta::GetSolutionFromInputsEnthalpy{{{1*/
    4213 void  Penta::GetSolutionFromInputsEnthalpy(Vec solution){
     4213void  Penta::GetSolutionFromInputsEnthalpy(Vector* solution){
    42144214
    42154215        const int    numdof=NDOF1*NUMVERTICES;
    … …  
    42344234
    42354235        /*Add value to global vector*/
    4236         VecSetValues(solution,numdof,doflist,(const double*)values,INSERT_VALUES);
     4236        solution->SetValues(numdof,doflist,values,INSERT_VALUES);
    42374237
    42384238        /*Free ressources:*/
    … …  
    78407840/*}}}*/
    78417841/*FUNCTION Penta::GetSolutionFromInputsDiagnosticHoriz{{{1*/
    7842 void  Penta::GetSolutionFromInputsDiagnosticHoriz(Vec solution){
     7842void  Penta::GetSolutionFromInputsDiagnosticHoriz(Vector* solution){
    78437843
    78447844        const int    numdof=NDOF2*NUMVERTICES;
    … …  
    78747874
    78757875        /*Add value to global vector*/
    7876         VecSetValues(solution,numdof,doflist,(const double*)values,INSERT_VALUES);
     7876        solution->SetValues(numdof,doflist,values,INSERT_VALUES);
    78777877
    78787878        /*Free ressources:*/
    … …  
    78827882/*}}}*/
    78837883/*FUNCTION Penta::GetSolutionFromInputsDiagnosticHutter{{{1*/
    7884 void  Penta::GetSolutionFromInputsDiagnosticHutter(Vec solution){
     7884void  Penta::GetSolutionFromInputsDiagnosticHutter(Vector* solution){
    78857885
    78867886        const int    numdof=NDOF2*NUMVERTICES;
    … …  
    79107910
    79117911        /*Add value to global vector*/
    7912         VecSetValues(solution,numdof,doflist,(const double*)values,INSERT_VALUES);
     7912        solution->SetValues(numdof,doflist,values,INSERT_VALUES);
    79137913
    79147914        /*Free ressources:*/
    … …  
    79187918/*}}}*/
    79197919/*FUNCTION Penta::GetSolutionFromInputsDiagnosticVert{{{1*/
    7920 void  Penta::GetSolutionFromInputsDiagnosticVert(Vec solution){
     7920void  Penta::GetSolutionFromInputsDiagnosticVert(Vector* solution){
    79217921
    79227922        const int    numdof=NDOF1*NUMVERTICES;
    … …  
    79437943
    79447944        /*Add value to global vector*/
    7945         VecSetValues(solution,numdof,doflist,(const double*)values,INSERT_VALUES);
     7945        solution->SetValues(numdof,doflist,values,INSERT_VALUES);
    79467946
    79477947        /*Free ressources:*/
    … …  
    79517951/*}}}*/
    79527952/*FUNCTION Penta::GetSolutionFromInputsDiagnosticStokes{{{1*/
    7953 void  Penta::GetSolutionFromInputsDiagnosticStokes(Vec solution){
     7953void  Penta::GetSolutionFromInputsDiagnosticStokes(Vector* solution){
    79547954
    79557955        const int    numdof=NDOF4*NUMVERTICES;
    … …  
    79887988
    79897989        /*Add value to global vector*/
    7990         VecSetValues(solution,numdof,doflist,(const double*)values,INSERT_VALUES);
     7990        solution->SetValues(numdof,doflist,values,INSERT_VALUES);
    79917991
    79927992        /*Free ressources:*/

  • issm/trunk-jpl/src/c/objects/Elements/Penta.h

    r11656 r11679  
    8686                void   Configure(Elements* elements,Loads* loads,DataSet* nodes,Materials* materials,Parameters* parameters);
    8787                void   SetCurrentConfiguration(Elements* elements,Loads* loads,DataSet* nodes,Materials* materials,Parameters* parameters);
    88                 void   CreateKMatrix(Mat Kff, Mat Kfs,Vec df);
    89                 void   CreatePVector(Vec pf);
    90                 void   CreateJacobianMatrix(Mat Jff);
     88                void   CreateKMatrix(Matrix* Kff, Matrix* Kfs,Vector* df);
     89                void   CreatePVector(Vector* pf);
     90                void   CreateJacobianMatrix(Matrix* Jff);
    9191                void   DeleteResults(void);
    9292                int    GetNodeIndex(Node* node);
    93                 void   GetSolutionFromInputs(Vec solution);
     93                void   GetSolutionFromInputs(Vector* solution);
    9494                double GetZcoord(GaussPenta* gauss);
    9595                void   GetVectorFromInputs(Vec vector,int name_enum);
    … …  
    183183                void    GetInputValue(double* pvalue,Node* node,int enumtype);
    184184                void      GetPhi(double* phi, double*  epsilon, double viscosity);
    185                 void      GetSolutionFromInputsEnthalpy(Vec solutiong);
     185                void      GetSolutionFromInputsEnthalpy(Vector* solutiong);
    186186                double  GetStabilizationParameter(double u, double v, double w, double diameter, double kappa);
    187187                void    GetStrainRate3dPattyn(double* epsilon,double* xyz_list, GaussPenta* gauss, Input* vx_input, Input* vy_input);
    … …  
    250250                void           InputUpdateFromSolutionDiagnosticVert( double* solutiong);
    251251                void           InputUpdateFromSolutionDiagnosticStokes( double* solutiong);
    252                 void             GetSolutionFromInputsDiagnosticHoriz(Vec solutiong);
    253                 void             GetSolutionFromInputsDiagnosticHutter(Vec solutiong);
    254                 void             GetSolutionFromInputsDiagnosticStokes(Vec solutiong);
    255                 void             GetSolutionFromInputsDiagnosticVert(Vec solutiong);
     252                void             GetSolutionFromInputsDiagnosticHoriz(Vector* solutiong);
     253                void             GetSolutionFromInputsDiagnosticHutter(Vector* solutiong);
     254                void             GetSolutionFromInputsDiagnosticStokes(Vector* solutiong);
     255                void             GetSolutionFromInputsDiagnosticVert(Vector* solutiong);
    256256                ElementVector* CreatePVectorCouplingMacAyealStokes(void);
    257257                ElementVector* CreatePVectorCouplingMacAyealStokesViscous(void);
    … …  
    307307                ElementVector* CreatePVectorThermalShelf(void);
    308308                ElementVector* CreatePVectorThermalSheet(void);
    309                 void           GetSolutionFromInputsThermal(Vec solutiong);
     309                void           GetSolutionFromInputsThermal(Vector* solutiong);
    310310                void           InputUpdateFromSolutionThermal( double* solutiong);
    311311                void           InputUpdateFromSolutionEnthalpy( double* solutiong);

  • issm/trunk-jpl/src/c/objects/Elements/Tria.cpp

    r11656 r11679  
    319319/*}}}*/
    320320/*FUNCTION Tria::CreateKMatrix {{{1*/
    321 void  Tria::CreateKMatrix(Mat Kff, Mat Kfs,Vec df){
     321void  Tria::CreateKMatrix(Matrix* Kff, Matrix* Kfs,Vector* df){
    322322
    323323        /*retreive parameters: */
    … …  
    672672/*}}}*/
    673673/*FUNCTION Tria::CreatePVector {{{1*/
    674 void  Tria::CreatePVector(Vec pf){
     674void  Tria::CreatePVector(Vector* pf){
    675675
    676676        /*retrive parameters: */
    … …  
    892892/*}}}*/
    893893/*FUNCTION Tria::CreateJacobianMatrix{{{1*/
    894 void  Tria::CreateJacobianMatrix(Mat Jff){
     894void  Tria::CreateJacobianMatrix(Matrix* Jff){
    895895
    896896        /*retrieve parameters: */
    … …  
    12741274/*}}}*/
    12751275/*FUNCTION Tria::GetSolutionFromInputs{{{1*/
    1276 void  Tria::GetSolutionFromInputs(Vec solution){
     1276void  Tria::GetSolutionFromInputs(Vector* solution){
    12771277
    12781278        /*retrive parameters: */
    … …  
    31393139/*}}}*/
    31403140/*FUNCTION Tria::GetSolutionFromInputsDiagnosticHoriz{{{1*/
    3141 void  Tria::GetSolutionFromInputsDiagnosticHoriz(Vec solution){
     3141void  Tria::GetSolutionFromInputsDiagnosticHoriz(Vector* solution){
    31423142
    31433143        const int    numdof=NDOF2*NUMVERTICES;
    … …  
    31703170        }
    31713171
    3172         VecSetValues(solution,numdof,doflist,(const double*)values,INSERT_VALUES);
     3172        solution->SetValues(numdof,doflist,values,INSERT_VALUES);
    31733173
    31743174        /*Free ressources:*/
    … …  
    31783178/*}}}*/
    31793179/*FUNCTION Tria::GetSolutionFromInputsDiagnosticHutter{{{1*/
    3180 void  Tria::GetSolutionFromInputsDiagnosticHutter(Vec solution){
     3180void  Tria::GetSolutionFromInputsDiagnosticHutter(Vector* solution){
    31813181
    31823182        const int    numdof=NDOF2*NUMVERTICES;
    … …  
    32093209        }
    32103210
    3211         VecSetValues(solution,numdof,doflist,(const double*)values,INSERT_VALUES);
     3211        solution->SetValues(numdof,doflist,values,INSERT_VALUES);
    32123212
    32133213        /*Free ressources:*/
    … …  
    48814881}
    48824882/*}}}*/
     4883
    48834884/*FUNCTION Tria::InputUpdateFromSolutionAdjointBalancethickness {{{1*/
    48844885void  Tria::InputUpdateFromSolutionAdjointBalancethickness(double* solution){
    … …  
    51855186/*}}}*/
    51865187/*FUNCTION Tria::GetSolutionFromInputsHydrology{{{1*/
    5187 void  Tria::GetSolutionFromInputsHydrology(Vec solution){
     5188void  Tria::GetSolutionFromInputsHydrology(Vector* solution){
    51885189
    51895190        const int    numdof=NDOF1*NUMVERTICES;
    … …  
    52135214        }
    52145215
    5215         VecSetValues(solution,numdof,doflist,(const double*)values,INSERT_VALUES);
     5216        solution->SetValues(numdof,doflist,values,INSERT_VALUES);
    52165217
    52175218        /*Free ressources:*/

  • issm/trunk-jpl/src/c/objects/Elements/Tria.h

    r11656 r11679  
    8282                void   Configure(Elements* elements,Loads* loads,DataSet* nodes,Materials* materials,Parameters* parameters);
    8383                void   SetCurrentConfiguration(Elements* elements,Loads* loads,DataSet* nodes,Materials* materials,Parameters* parameters);
    84                 void   CreateKMatrix(Mat Kff, Mat Kfs,Vec df);
    85                 void   CreatePVector(Vec pf);
    86                 void   CreateJacobianMatrix(Mat Jff);
     84                void   CreateKMatrix(Matrix* Kff, Matrix* Kfs,Vector* df);
     85                void   CreatePVector(Vector* pf);
     86                void   CreateJacobianMatrix(Matrix* Jff);
    8787                int    GetNodeIndex(Node* node);
    8888                int    Sid();
    … …  
    9292                bool   IsNodeOnShelfFromFlags(double* flags);
    9393                bool   IsOnWater();
    94                 void   GetSolutionFromInputs(Vec solution);
     94                void   GetSolutionFromInputs(Vector* solution);
    9595                void   GetVectorFromInputs(Vec vector, int name_enum);
    9696                void   GetVectorFromResults(Vec vector,int offset,int interp);
    … …  
    210210                ElementVector* CreatePVectorDiagnosticHutter(void);
    211211                ElementMatrix* CreateJacobianDiagnosticMacayeal(void);
    212                 void      GetSolutionFromInputsDiagnosticHoriz(Vec solution);
    213                 void      GetSolutionFromInputsDiagnosticHutter(Vec solution);
     212                void      GetSolutionFromInputsDiagnosticHoriz(Vector* solution);
     213                void      GetSolutionFromInputsDiagnosticHutter(Vector* solution);
    214214                void      InputUpdateFromSolutionDiagnosticHoriz( double* solution);
    215215                void      InputUpdateFromSolutionDiagnosticHutter( double* solution);
    … …  
    230230                ElementVector* CreatePVectorHydrology(void);
    231231                void      CreateHydrologyWaterVelocityInput(void);
    232                 void      GetSolutionFromInputsHydrology(Vec solution);
     232                void      GetSolutionFromInputsHydrology(Vector* solution);
    233233                void      InputUpdateFromSolutionHydrology(double* solution);
    234234                #endif

  • issm/trunk-jpl/src/c/objects/Loads/Icefront.cpp

    r11332 r11679  
    316316/*}}}*/
    317317/*FUNCTION Icefront::CreateKMatrix {{{1*/
    318 void  Icefront::CreateKMatrix(Mat Kff, Mat Kfs){
     318void  Icefront::CreateKMatrix(Matrix* Kff, Matrix* Kfs){
    319319
    320320        /*No stiffness loads applied, do nothing: */
    … …  
    324324/*}}}*/
    325325/*FUNCTION Icefront::CreatePVector {{{1*/
    326 void  Icefront::CreatePVector(Vec pf){
     326void  Icefront::CreatePVector(Vector* pf){
    327327
    328328        /*Checks in debugging mode*/
    … …  
    362362/*}}}*/
    363363/*FUNCTION Icefront::CreateJacobianMatrix{{{1*/
    364 void  Icefront::CreateJacobianMatrix(Mat Jff){
     364void  Icefront::CreateJacobianMatrix(Matrix* Jff){
    365365        this->CreateKMatrix(Jff,NULL);
    366366}
    367367/*}}}1*/
    368368/*FUNCTION Icefront::PenaltyCreateKMatrix {{{1*/
    369 void  Icefront::PenaltyCreateKMatrix(Mat Kff, Mat Kfs, double kmax){
     369void  Icefront::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs, double kmax){
    370370        /*do nothing: */
    371371        return;
    … …  
    373373/*}}}*/
    374374/*FUNCTION Icefront::PenaltyCreatePVector{{{1*/
    375 void  Icefront::PenaltyCreatePVector(Vec pf,double kmax){
     375void  Icefront::PenaltyCreatePVector(Vector* pf,double kmax){
    376376        /*do nothing: */
    377377        return;
    … …  
    379379/*}}}*/
    380380/*FUNCTION Icefront::PenaltyCreateJacobianMatrix{{{1*/
    381 void  Icefront::PenaltyCreateJacobianMatrix(Mat Jff,double kmax){
     381void  Icefront::PenaltyCreateJacobianMatrix(Matrix* Jff,double kmax){
    382382        this->PenaltyCreateKMatrix(Jff,NULL,kmax);
    383383}

  • issm/trunk-jpl/src/c/objects/Loads/Icefront.h

    r11332 r11679  
    7474                void  Configure(Elements* elements,Loads* loads,Nodes* nodes,Vertices* vertices,Materials* materials,Parameters* parameters);
    7575                void  SetCurrentConfiguration(Elements* elements,Loads* loads,Nodes* nodes,Vertices* vertices,Materials* materials,Parameters* parameters);
    76                 void  CreateKMatrix(Mat Kff, Mat Kfs);
    77                 void  CreatePVector(Vec pf);
    78                 void  CreateJacobianMatrix(Mat Jff);
    79                 void  PenaltyCreateKMatrix(Mat Kff, Mat kfs, double kmax);
    80                 void  PenaltyCreatePVector(Vec pf, double kmax);
    81                 void  PenaltyCreateJacobianMatrix(Mat Jff,double kmax);
     76                void  CreateKMatrix(Matrix* Kff, Matrix* Kfs);
     77                void  CreatePVector(Vector* pf);
     78                void  CreateJacobianMatrix(Matrix* Jff);
     79                void  PenaltyCreateKMatrix(Matrix* Kff, Matrix* kfs, double kmax);
     80                void  PenaltyCreatePVector(Vector* pf, double kmax);
     81                void  PenaltyCreateJacobianMatrix(Matrix* Jff,double kmax);
    8282                bool  InAnalysis(int analysis_type);
    8383                /*}}}*/

  • issm/trunk-jpl/src/c/objects/Loads/Load.h

    r11327 r11679  
    1212/*{{{1*/
    1313class Object;
     14class Matrix;
     15class Vector;
    1416
    1517#include "../Object.h"
    … …  
    2729                virtual void  Configure(Elements* elements,Loads* loads,Nodes* nodes,Vertices* vertices,Materials* materials,Parameters* parameters)=0;
    2830                virtual void  SetCurrentConfiguration(Elements* elements,Loads* loads,Nodes* nodes,Vertices* vertices,Materials* materials,Parameters* parameters)=0;
    29                 virtual void  CreateKMatrix(Mat Kff, Mat Kfs)=0;
    30                 virtual void  CreatePVector(Vec pf)=0;
    31                 virtual void  CreateJacobianMatrix(Mat Jff)=0;
    32                 virtual void  PenaltyCreateJacobianMatrix(Mat Jff,double kmax)=0;
    33                 virtual void  PenaltyCreateKMatrix(Mat Kff, Mat Kfs, double kmax)=0;
    34                 virtual void  PenaltyCreatePVector(Vec pf, double kmax)=0;
     31                virtual void  CreateKMatrix(Matrix* Kff, Matrix* Kfs)=0;
     32                virtual void  CreatePVector(Vector* pf)=0;
     33                virtual void  CreateJacobianMatrix(Matrix* Jff)=0;
     34                virtual void  PenaltyCreateJacobianMatrix(Matrix* Jff,double kmax)=0;
     35                virtual void  PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs, double kmax)=0;
     36                virtual void  PenaltyCreatePVector(Vector* pf, double kmax)=0;
    3537                virtual bool  InAnalysis(int analysis_type)=0;
    3638                /*}}}*/

  • issm/trunk-jpl/src/c/objects/Loads/Numericalflux.cpp

    r10135 r11679  
    334334/*}}}*/
    335335/*FUNCTION Numericalflux::CreateKMatrix {{{1*/
    336 void  Numericalflux::CreateKMatrix(Mat Kff, Mat Kfs){
     336void  Numericalflux::CreateKMatrix(Matrix* Kff, Matrix* Kfs){
    337337
    338338        /*recover some parameters*/
    … …  
    365365/*}}}*/
    366366/*FUNCTION Numericalflux::CreatePVector {{{1*/
    367 void  Numericalflux::CreatePVector(Vec pf){
     367void  Numericalflux::CreatePVector(Vector* pf){
    368368
    369369        /*recover some parameters*/
    … …  
    395395/*}}}*/
    396396/*FUNCTION Numericalflux::PenaltyCreateKMatrix {{{1*/
    397 void  Numericalflux::PenaltyCreateKMatrix(Mat Kff, Mat Kfs,double kmax){
     397void  Numericalflux::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs,double kmax){
    398398
    399399        /*No stiffness loads applied, do nothing: */
    … …  
    403403/*}}}*/
    404404/*FUNCTION Numericalflux::PenaltyCreatePVector{{{1*/
    405 void  Numericalflux::PenaltyCreatePVector(Vec pf,double kmax){
     405void  Numericalflux::PenaltyCreatePVector(Vector* pf,double kmax){
    406406
    407407        /*No penalty loads applied, do nothing: */

  • issm/trunk-jpl/src/c/objects/Loads/Numericalflux.h

    r11327 r11679  
    7070                void  Configure(Elements* elements,Loads* loads,Nodes* nodes,Vertices* vertices,Materials* materials,Parameters* parameters);
    7171                void  SetCurrentConfiguration(Elements* elements,Loads* loads,Nodes* nodes,Vertices* vertices,Materials* materials,Parameters* parameters);
    72                 void  CreateKMatrix(Mat Kff, Mat Kfs);
    73                 void  CreatePVector(Vec pf);
    74                 void  CreateJacobianMatrix(Mat Jff){_error_("Not implemented yet");};
    75                 void  PenaltyCreateJacobianMatrix(Mat Jff,double kmax){_error_("Not implemented yet");};
    76                 void  PenaltyCreateKMatrix(Mat Kff, Mat kfs, double kmax);
    77                 void  PenaltyCreatePVector(Vec pf, double kmax);
     72                void  CreateKMatrix(Matrix* Kff, Matrix* Kfs);
     73                void  CreatePVector(Vector* pf);
     74                void  CreateJacobianMatrix(Matrix* Jff){_error_("Not implemented yet");};
     75                void  PenaltyCreateJacobianMatrix(Matrix* Jff,double kmax){_error_("Not implemented yet");};
     76                void  PenaltyCreateKMatrix(Matrix* Kff, Matrix* kfs, double kmax);
     77                void  PenaltyCreatePVector(Vector* pf, double kmax);
    7878                bool  InAnalysis(int analysis_type);
    7979                /*}}}*/

  • issm/trunk-jpl/src/c/objects/Loads/Pengrid.cpp

    r11656 r11679  
    295295/*}}}1*/
    296296/*FUNCTION Pengrid::CreateKMatrix {{{1*/
    297 void  Pengrid::CreateKMatrix(Mat Kff, Mat Kfs){
     297void  Pengrid::CreateKMatrix(Matrix* Kff, Matrix* Kfs){
    298298
    299299        /*No loads applied, do nothing: */
    … …  
    303303/*}}}1*/
    304304/*FUNCTION Pengrid::CreatePVector {{{1*/
    305 void  Pengrid::CreatePVector(Vec pf){
     305void  Pengrid::CreatePVector(Vector* pf){
    306306
    307307        /*No loads applied, do nothing: */
    … …  
    311311/*}}}1*/
    312312/*FUNCTION Pengrid::PenaltyCreateMatrix {{{1*/
    313 void  Pengrid::PenaltyCreateKMatrix(Mat Kff, Mat Kfs,double kmax){
     313void  Pengrid::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs,double kmax){
    314314
    315315        /*Retrieve parameters: */
    … …  
    344344/*}}}1*/
    345345/*FUNCTION Pengrid::PenaltyCreatePVector {{{1*/
    346 void  Pengrid::PenaltyCreatePVector(Vec pf,double kmax){
     346void  Pengrid::PenaltyCreatePVector(Vector* pf,double kmax){
    347347
    348348        /*Retrieve parameters: */

  • issm/trunk-jpl/src/c/objects/Loads/Pengrid.h

    r11656 r11679  
    7575                void  Configure(Elements* elements,Loads* loads,Nodes* nodes,Vertices* vertices,Materials* materials,Parameters* parameters);
    7676                void  SetCurrentConfiguration(Elements* elements,Loads* loads,Nodes* nodes,Vertices* vertices,Materials* materials,Parameters* parameters);
    77                 void  CreateKMatrix(Mat Kff, Mat Kfs);
    78                 void  CreatePVector(Vec pf);
    79                 void  CreateJacobianMatrix(Mat Jff){_error_("Not implemented yet");};
    80                 void  PenaltyCreateJacobianMatrix(Mat Jff,double kmax){_error_("Not implemented yet");};
    81                 void  PenaltyCreateKMatrix(Mat Kff, Mat kfs, double kmax);
    82                 void  PenaltyCreatePVector(Vec pf, double kmax);
     77                void  CreateKMatrix(Matrix* Kff, Matrix* Kfs);
     78                void  CreatePVector(Vector* pf);
     79                void  CreateJacobianMatrix(Matrix* Jff){_error_("Not implemented yet");};
     80                void  PenaltyCreateJacobianMatrix(Matrix* Jff,double kmax){_error_("Not implemented yet");};
     81                void  PenaltyCreateKMatrix(Matrix* Kff, Matrix* kfs, double kmax);
     82                void  PenaltyCreatePVector(Vector* pf, double kmax);
    8383                bool  InAnalysis(int analysis_type);
    8484                /*}}}*/

  • issm/trunk-jpl/src/c/objects/Loads/Penpair.cpp

    r11332 r11679  
    199199/*}}}1*/
    200200/*FUNCTION Penpair::CreateKMatrix {{{1*/
    201 void  Penpair::CreateKMatrix(Mat Kff, Mat Kfs){
     201void  Penpair::CreateKMatrix(Matrix* Kff, Matrix* Kfs){
    202202        /*If you code this piece, don't forget that a penalty will be inactive if it is dealing with clone nodes*/
    203203        /*No loads applied, do nothing: */
    … …  
    207207/*}}}1*/
    208208/*FUNCTION Penpair::CreatePVector {{{1*/
    209 void  Penpair::CreatePVector(Vec pf){
     209void  Penpair::CreatePVector(Vector* pf){
    210210
    211211        /*No loads applied, do nothing: */
    … …  
    215215/*}}}1*/
    216216/*FUNCTION Penpair::CreateJacobianMatrix{{{1*/
    217 void  Penpair::CreateJacobianMatrix(Mat Jff){
     217void  Penpair::CreateJacobianMatrix(Matrix* Jff){
    218218        this->CreateKMatrix(Jff,NULL);
    219219}
    220220/*}}}1*/
    221221/*FUNCTION Penpair::PenaltyCreateKMatrix {{{1*/
    222 void  Penpair::PenaltyCreateKMatrix(Mat Kff, Mat Kfs,double kmax){
     222void  Penpair::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs,double kmax){
    223223
    224224        /*Retrieve parameters: */
    … …  
    246246/*}}}1*/
    247247/*FUNCTION Penpair::PenaltyCreatePVector {{{1*/
    248 void  Penpair::PenaltyCreatePVector(Vec pf,double kmax){
     248void  Penpair::PenaltyCreatePVector(Vector* pf,double kmax){
    249249        /*No loads applied, do nothing: */
    250250        return;
    … …  
    252252/*}}}1*/
    253253/*FUNCTION Penpair::PenaltyCreateJacobianMatrix{{{1*/
    254 void  Penpair::PenaltyCreateJacobianMatrix(Mat Jff,double kmax){
     254void  Penpair::PenaltyCreateJacobianMatrix(Matrix* Jff,double kmax){
    255255        this->PenaltyCreateKMatrix(Jff,NULL,kmax);
    256256}

  • issm/trunk-jpl/src/c/objects/Loads/Penpair.h

    r11327 r11679  
    6262                void  Configure(Elements* elements,Loads* loads,Nodes* nodes,Vertices* vertices,Materials* materials,Parameters* parameters);
    6363                void  SetCurrentConfiguration(Elements* elements,Loads* loads,Nodes* nodes,Vertices* vertices,Materials* materials,Parameters* parameters);
    64                 void  CreateKMatrix(Mat Kff, Mat Kfs);
    65                 void  CreatePVector(Vec pf);
    66                 void  CreateJacobianMatrix(Mat Jff);
    67                 void  PenaltyCreateKMatrix(Mat Kff,Mat Kfs,double kmax);
    68                 void  PenaltyCreatePVector(Vec pf, double kmax);
    69                 void  PenaltyCreateJacobianMatrix(Mat Jff,double kmax);
     64                void  CreateKMatrix(Matrix* Kff, Matrix* Kfs);
     65                void  CreatePVector(Vector* pf);
     66                void  CreateJacobianMatrix(Matrix* Jff);
     67                void  PenaltyCreateKMatrix(Matrix* Kff,Matrix* Kfs,double kmax);
     68                void  PenaltyCreatePVector(Vector* pf, double kmax);
     69                void  PenaltyCreateJacobianMatrix(Matrix* Jff,double kmax);
    7070                bool  InAnalysis(int analysis_type);
    7171                /*}}}*/

  • issm/trunk-jpl/src/c/objects/Loads/Riftfront.cpp

    r11294 r11679  
    428428/*}}}*/
    429429/*FUNCTION Riftfront::PenaltyCreateKMatrix {{{1*/
    430 void  Riftfront::PenaltyCreateKMatrix(Mat Kff, Mat Kfs,double kmax){
     430void  Riftfront::PenaltyCreateKMatrix(Matrix* Kff, Matrix* Kfs,double kmax){
    431431
    432432        /*Retrieve parameters: */
    … …  
    454454/*}}}1*/
    455455/*FUNCTION Riftfront::PenaltyCreatePVector {{{1*/
    456 void  Riftfront::PenaltyCreatePVector(Vec pf,double kmax){
     456void  Riftfront::PenaltyCreatePVector(Vector* pf,double kmax){
    457457
    458458        /*Retrieve parameters: */
    … …  
    480480/*}}}1*/
    481481/*FUNCTION Riftfront::CreateKMatrix {{{1*/
    482 void  Riftfront::CreateKMatrix(Mat Kff, Mat Kfs){
     482void  Riftfront::CreateKMatrix(Matrix* Kff, Matrix* Kfs){
    483483        /*do nothing: */
    484484        return;
    … …  
    486486/*}}}1*/
    487487/*FUNCTION Riftfront::CreatePVector {{{1*/
    488 void  Riftfront::CreatePVector(Vec pf){
     488void  Riftfront::CreatePVector(Vector* pf){
    489489        /*do nothing: */
    490490        return;

  • issm/trunk-jpl/src/c/objects/Loads/Riftfront.h

    r11327 r11679  
    8282                void  Configure(Elements* elements,Loads* loads,Nodes* nodes,Vertices* vertices,Materials* materials,Parameters* parameters);
    8383                void  SetCurrentConfiguration(Elements* elements,Loads* loads,Nodes* nodes,Vertices* vertices,Materials* materials,Parameters* parameters);
    84                 void  CreateKMatrix(Mat Kff, Mat Kfs);
    85                 void  CreatePVector(Vec pf);
    86                 void  CreateJacobianMatrix(Mat Jff){_error_("Not implemented yet");};
    87                 void  PenaltyCreateJacobianMatrix(Mat Jff,double kmax){_error_("Not implemented yet");};
    88                 void  PenaltyCreateKMatrix(Mat Kff, Mat kfs, double kmax);
    89                 void  PenaltyCreatePVector(Vec pf, double kmax);
     84                void  CreateKMatrix(Matrix* Kff, Matrix* Kfs);
     85                void  CreatePVector(Vector* pf);
     86                void  CreateJacobianMatrix(Matrix* Jff){_error_("Not implemented yet");};
     87                void  PenaltyCreateJacobianMatrix(Matrix* Jff,double kmax){_error_("Not implemented yet");};
     88                void  PenaltyCreateKMatrix(Matrix* Kff, Matrix* kfs, double kmax);
     89                void  PenaltyCreatePVector(Vector* pf, double kmax);
    9090                bool  InAnalysis(int analysis_type);
    9191                /*}}}*/

  • issm/trunk-jpl/src/c/objects/Node.cpp

    r11001 r11679  
    594594/*}}}*/
    595595/*FUNCTION Node::CreateNodalConstraints{{{1*/
    596 void  Node::CreateNodalConstraints(Vec ys){
     596void  Node::CreateNodalConstraints(Vector* ys){
    597597
    598598        int i;
    … …  
    613613               
    614614                /*Add values into constraint vector: */
    615                 VecSetValues(ys,this->indexing.ssize,this->indexing.sdoflist,values,INSERT_VALUES);
     615                ys->SetValues(this->indexing.ssize,this->indexing.sdoflist,values,INSERT_VALUES);
    616616        }
    617617
    … …  
    875875/*}}}*/
    876876/*FUNCTION Node::VecMerge {{{1*/
    877 void   Node::VecMerge(Vec ug, double* vector_serial,int setenum){
     877void   Node::VecMerge(Vector* ug, double* vector_serial,int setenum){
    878878
    879879        double* values=NULL;
    … …  
    897897
    898898                        /*Add values into ug: */
    899                         VecSetValues(ug,this->indexing.fsize,indices,(const double*)values,INSERT_VALUES);
     899                        ug->SetValues(this->indexing.fsize,indices,values,INSERT_VALUES);
    900900                }
    901901        }
    … …  
    915915
    916916                        /*Add values into ug: */
    917                         VecSetValues(ug,this->indexing.ssize,indices,(const double*)values,INSERT_VALUES);
     917                        ug->SetValues(this->indexing.ssize,indices,values,INSERT_VALUES);
    918918                }
    919919        }
    … …  
    926926/*}}}*/
    927927/*FUNCTION Node::VecReduce {{{1*/
    928 void   Node::VecReduce(Vec vector, double* ug_serial,int setenum){
     928void   Node::VecReduce(Vector* vector, double* ug_serial,int setenum){
    929929
    930930        double* values=NULL;
    … …  
    945945
    946946                        /*Add values into ug: */
    947                         VecSetValues(vector,this->indexing.fsize,this->indexing.fdoflist,(const double*)values,INSERT_VALUES);
     947                        vector->SetValues(this->indexing.fsize,this->indexing.fdoflist,values,INSERT_VALUES);
    948948                }
    949949        }
    … …  
    961961
    962962                        /*Add values into ug: */
    963                         VecSetValues(vector,this->indexing.ssize,this->indexing.sdoflist,(const double*)values,INSERT_VALUES);
     963                        vector->SetValues(this->indexing.ssize,this->indexing.sdoflist,values,INSERT_VALUES);
    964964                }
    965965        }

  • issm/trunk-jpl/src/c/objects/Node.h

    r10576 r11679  
    1616class  DataSet;
    1717class  Vertices;
     18class  Vector;
     19class  Matrix;
    1820#include "./Update.h"
    1921/*}}}*/
    … …  
    6769                /*Node numerical routines {{{1*/
    6870                void   Configure(DataSet* nodes,Vertices* vertices);
    69                 void   CreateNodalConstraints(Vec ys);
     71                void   CreateNodalConstraints(Vector* ys);
    7072                void   SetCurrentConfiguration(DataSet* nodes,Vertices* vertices);
    7173                int    Sid(void);
    … …  
    101103                int    IsGrounded();
    102104                void   UpdateSpcs(double* ys);
    103                 void   VecMerge(Vec ug, double* vector_serial,int setnum);
    104                 void   VecReduce(Vec vector, double* ug_serial,int setnum);
     105                void   VecMerge(Vector* ug, double* vector_serial,int setenum);
     106                void   VecReduce(Vector* vector, double* ug_serial,int setnum);
    105107               
    106108                /*}}}*/

  • issm/trunk-jpl/src/c/objects/Numerics/ElementMatrix.cpp

    r11327 r11679  
    245245
    246246/*ElementMatrix specific routines: */
    247 /*FUNCTION ElementMatrix::AddToGlobal(Mat Kff, Mat Kfs){{{1*/
    248 void ElementMatrix::AddToGlobal(Mat Kff, Mat Kfs){
     247/*FUNCTION ElementMatrix::AddToGlobal(Matrix* Kff, Matrix* Kfs){{{1*/
     248void ElementMatrix::AddToGlobal(Matrix* Kff, Matrix* Kfs){
    249249
    250250        int i,j;
    … …  
    260260        this->CheckConsistency();
    261261
    262         if(this->dofsymmetrical){
    263                 /*only use row dofs to add values into global matrices: */
    264                
    265                 if(this->row_fsize){
    266                         /*first, retrieve values that are in the f-set from the g-set values matrix: */
    267                         localvalues=(double*)xmalloc(this->row_fsize*this->row_fsize*sizeof(double));
    268                         for(i=0;i<this->row_fsize;i++){
    269                                 for(j=0;j<this->row_fsize;j++){
    270                                         *(localvalues+this->row_fsize*i+j)=*(this->values+this->ncols*this->row_flocaldoflist[i]+this->row_flocaldoflist[j]);
     262        #ifdef _HAVE_PETSC_
     263                if(this->dofsymmetrical){
     264                        /*only use row dofs to add values into global matrices: */
     265                       
     266                        if(this->row_fsize){
     267                                /*first, retrieve values that are in the f-set from the g-set values matrix: */
     268                                localvalues=(double*)xmalloc(this->row_fsize*this->row_fsize*sizeof(double));
     269                                for(i=0;i<this->row_fsize;i++){
     270                                        for(j=0;j<this->row_fsize;j++){
     271                                                *(localvalues+this->row_fsize*i+j)=*(this->values+this->ncols*this->row_flocaldoflist[i]+this->row_flocaldoflist[j]);
     272                                        }
    271273                                }
     274                                /*add local values into global  matrix, using the fglobaldoflist: */
     275                                MatSetValues(Kff->matrix,this->row_fsize,this->row_fglobaldoflist,this->row_fsize,this->row_fglobaldoflist,(const double*)localvalues,ADD_VALUES);
     276
     277                                /*Free ressources:*/
     278                                xfree((void**)&localvalues);
    272279                        }
    273                         /*add local values into global  matrix, using the fglobaldoflist: */
    274                         MatSetValues(Kff,this->row_fsize,this->row_fglobaldoflist,this->row_fsize,this->row_fglobaldoflist,(const double*)localvalues,ADD_VALUES);
    275 
    276                         /*Free ressources:*/
    277                         xfree((void**)&localvalues);
    278                 }
    279 
    280 
    281                 if((this->row_ssize!=0) && (this->row_fsize!=0)){
    282                         /*first, retrieve values that are in the f and s-set from the g-set values matrix: */
    283                         localvalues=(double*)xmalloc(this->row_fsize*this->row_ssize*sizeof(double));
    284                         for(i=0;i<this->row_fsize;i++){
    285                                 for(j=0;j<this->row_ssize;j++){
    286                                         *(localvalues+this->row_ssize*i+j)=*(this->values+this->ncols*this->row_flocaldoflist[i]+this->row_slocaldoflist[j]);
     280
     281
     282                        if((this->row_ssize!=0) && (this->row_fsize!=0)){
     283                                /*first, retrieve values that are in the f and s-set from the g-set values matrix: */
     284                                localvalues=(double*)xmalloc(this->row_fsize*this->row_ssize*sizeof(double));
     285                                for(i=0;i<this->row_fsize;i++){
     286                                        for(j=0;j<this->row_ssize;j++){
     287                                                *(localvalues+this->row_ssize*i+j)=*(this->values+this->ncols*this->row_flocaldoflist[i]+this->row_slocaldoflist[j]);
     288                                        }
    287289                                }
     290                                /*add local values into global  matrix, using the fglobaldoflist: */
     291                                MatSetValues(Kfs->matrix,this->row_fsize,this->row_fglobaldoflist,this->row_ssize,this->row_sglobaldoflist,(const double*)localvalues,ADD_VALUES);
     292
     293                                /*Free ressources:*/
     294                                xfree((void**)&localvalues);
    288295                        }
    289                         /*add local values into global  matrix, using the fglobaldoflist: */
    290                         MatSetValues(Kfs,this->row_fsize,this->row_fglobaldoflist,this->row_ssize,this->row_sglobaldoflist,(const double*)localvalues,ADD_VALUES);
    291 
    292                         /*Free ressources:*/
    293                         xfree((void**)&localvalues);
    294                 }
    295         }
    296         else{
    297                 _error_(" non dofsymmetrical matrix AddToGlobal routine not support yet!");
    298         }
    299 
    300 }
    301 /*}}}*/
    302 /*FUNCTION ElementMatrix::AddToGlobal(Mat Jff){{{1*/
    303 void ElementMatrix::AddToGlobal(Mat Jff){
     296                }
     297                else{
     298                        _error_(" non dofsymmetrical matrix AddToGlobal routine not support yet!");
     299                }
     300        #else
     301                _error_("not supported yet!");
     302        #endif
     303
     304}
     305/*}}}*/
     306/*FUNCTION ElementMatrix::AddToGlobal(Matrix* Jff){{{1*/
     307void ElementMatrix::AddToGlobal(Matrix* Jff){
    304308
    305309        int i,j;
    … …  
    312316        this->CheckConsistency();
    313317
    314         if(this->dofsymmetrical){
    315                 /*only use row dofs to add values into global matrices: */
    316 
    317                 if(this->row_fsize){
    318                         /*first, retrieve values that are in the f-set from the g-set values matrix: */
    319                         localvalues=(double*)xmalloc(this->row_fsize*this->row_fsize*sizeof(double));
    320                         for(i=0;i<this->row_fsize;i++){
    321                                 for(j=0;j<this->row_fsize;j++){
    322                                         *(localvalues+this->row_fsize*i+j)=*(this->values+this->ncols*this->row_flocaldoflist[i]+this->row_flocaldoflist[j]);
     318        #ifdef _HAVE_PETSC_
     319                if(this->dofsymmetrical){
     320                        /*only use row dofs to add values into global matrices: */
     321
     322                        if(this->row_fsize){
     323                                /*first, retrieve values that are in the f-set from the g-set values matrix: */
     324                                localvalues=(double*)xmalloc(this->row_fsize*this->row_fsize*sizeof(double));
     325                                for(i=0;i<this->row_fsize;i++){
     326                                        for(j=0;j<this->row_fsize;j++){
     327                                                *(localvalues+this->row_fsize*i+j)=*(this->values+this->ncols*this->row_flocaldoflist[i]+this->row_flocaldoflist[j]);
     328                                        }
    323329                                }
     330                                /*add local values into global  matrix, using the fglobaldoflist: */
     331                                MatSetValues(Jff->matrix,this->row_fsize,this->row_fglobaldoflist,this->row_fsize,this->row_fglobaldoflist,(const double*)localvalues,ADD_VALUES);
     332
     333                                /*Free ressources:*/
     334                                xfree((void**)&localvalues);
    324335                        }
    325                         /*add local values into global  matrix, using the fglobaldoflist: */
    326                         MatSetValues(Jff,this->row_fsize,this->row_fglobaldoflist,this->row_fsize,this->row_fglobaldoflist,(const double*)localvalues,ADD_VALUES);
    327 
    328                         /*Free ressources:*/
    329                         xfree((void**)&localvalues);
    330                 }
    331 
    332         }
    333         else{
    334                 _error_(" non dofsymmetrical matrix AddToGlobal routine not support yet!");
    335         }
     336
     337                }
     338                else{
     339                        _error_(" non dofsymmetrical matrix AddToGlobal routine not support yet!");
     340                }
     341        #else
     342                _error_("not supported yet!");
     343        #endif
    336344
    337345}

  • issm/trunk-jpl/src/c/objects/Numerics/ElementMatrix.h

    r11322 r11679  
    5858                /*}}}*/
    5959                /*ElementMatrix specific routines {{{1*/
    60                 void AddToGlobal(Mat Kff, Mat Kfs);
    61                 void AddToGlobal(Mat Jff);
     60                void AddToGlobal(Matrix* Kff, Matrix* Kfs);
     61                void AddToGlobal(Matrix* Jff);
    6262                void Echo(void);
    6363                void CheckConsistency(void);

  • issm/trunk-jpl/src/c/objects/Numerics/ElementVector.cpp

    r9320 r11679  
    160160
    161161/*ElementVector specific routines: */
    162 /*FUNCTION ElementVector::AddToGlobal(Vec pf){{{1*/
    163 void ElementVector::AddToGlobal(Vec pf){
     162/*FUNCTION ElementVector::AddToGlobal(Vector* pf){{{1*/
     163void ElementVector::AddToGlobal(Vector* pf){
    164164
    165165        int i;
    166166        double* localvalues=NULL;
    167167
    168         if(this->fsize){
    169                 /*first, retrieve values that are in the f-set from the g-set values vector: */
    170                 localvalues=(double*)xmalloc(this->fsize*sizeof(double));
    171                 for(i=0;i<this->fsize;i++){
    172                         localvalues[i]=this->values[this->flocaldoflist[i]];
    173                 }
    174                 /*add local values into global  vector, using the fglobaldoflist: */
    175                 VecSetValues(pf,this->fsize,this->fglobaldoflist,(const double*)localvalues,ADD_VALUES);
    176 
    177                 /*Free ressources:*/
    178                 xfree((void**)&localvalues);
    179         }
    180 }
    181 /*}}}*/
    182 /*FUNCTION ElementVector::InsertIntoGlobal(Vec pf){{{1*/
    183 void ElementVector::InsertIntoGlobal(Vec pf){
     168        #ifdef _HAVE_PETSC_
     169                if(this->fsize){
     170                        /*first, retrieve values that are in the f-set from the g-set values vector: */
     171                        localvalues=(double*)xmalloc(this->fsize*sizeof(double));
     172                        for(i=0;i<this->fsize;i++){
     173                                localvalues[i]=this->values[this->flocaldoflist[i]];
     174                        }
     175                        /*add local values into global  vector, using the fglobaldoflist: */
     176                        VecSetValues(pf->vector,this->fsize,this->fglobaldoflist,(const double*)localvalues,ADD_VALUES);
     177
     178                        /*Free ressources:*/
     179                        xfree((void**)&localvalues);
     180                }
     181        #else
     182                _error_("not supported yet!");
     183        #endif
     184       
     185}
     186/*}}}*/
     187/*FUNCTION ElementVector::InsertIntoGlobal(Vector* pf){{{1*/
     188void ElementVector::InsertIntoGlobal(Vector* pf){
    184189
    185190        int i;
    186191        double* localvalues=NULL;
    187192
    188         if(this->fsize){
    189                 /*first, retrieve values that are in the f-set from the g-set values vector: */
    190                 localvalues=(double*)xmalloc(this->fsize*sizeof(double));
    191                 for(i=0;i<this->fsize;i++){
    192                         localvalues[i]=this->values[this->flocaldoflist[i]];
    193                 }
    194                 /*add local values into global  vector, using the fglobaldoflist: */
    195                 VecSetValues(pf,this->fsize,this->fglobaldoflist,(const double*)localvalues,INSERT_VALUES);
    196 
    197                 /*Free ressources:*/
    198                 xfree((void**)&localvalues);
    199         }
     193        #ifdef _HAVE_PETSC_
     194                if(this->fsize){
     195                        /*first, retrieve values that are in the f-set from the g-set values vector: */
     196                        localvalues=(double*)xmalloc(this->fsize*sizeof(double));
     197                        for(i=0;i<this->fsize;i++){
     198                                localvalues[i]=this->values[this->flocaldoflist[i]];
     199                        }
     200                        /*add local values into global  vector, using the fglobaldoflist: */
     201                        VecSetValues(pf->vector,this->fsize,this->fglobaldoflist,(const double*)localvalues,INSERT_VALUES);
     202
     203                        /*Free ressources:*/
     204                        xfree((void**)&localvalues);
     205                }
     206        #else
     207                _error_("not supported yet!");
     208        #endif
     209
    200210}
    201211/*}}}*/

  • issm/trunk-jpl/src/c/objects/Numerics/ElementVector.h

    r8800 r11679  
    4040                /*}}}*/
    4141                /*ElementVector specific routines {{{1*/
    42                 void AddToGlobal(Vec pf);
    43                 void InsertIntoGlobal(Vec pf);
     42                void AddToGlobal(Vector* pf);
     43                void InsertIntoGlobal(Vector* pf);
    4444                void Echo(void);
    4545                void Init(ElementVector* pe);

  • issm/trunk-jpl/src/c/objects/Numerics/Matrix.cpp

    r11670 r11679  
    140140        _error_("not implemented yet!");
    141141        #endif
     142        return dataref;
    142143
    143144}
    144145/*}}}*/
    145146#endif
     147/*FUNCTION Matrix::Assemble{{{1*/
     148void Matrix::Assemble(void){
     149        #ifdef _HAVE_PETSC_
     150                MatAssemblyBegin(this->matrix,MAT_FINAL_ASSEMBLY);
     151                MatAssemblyEnd(this->matrix,MAT_FINAL_ASSEMBLY);
     152                #if _PETSC_MAJOR_ == 2
     153                        MatCompress(this->matrix);
     154                #endif
     155        #else
     156                /*do nothing:*/
     157        #endif
     158
     159}
     160/*}}}*/
     161/*FUNCTION Matrix::Norm{{{1*/
     162double Matrix::Norm(int norm_type){
     163       
     164        double norm=0;
     165        #ifdef _HAVE_PETSC_
     166            MatNorm(this->matrix,(NormType)norm_type,&norm);
     167        #else
     168                _error_("not implemented yet!");
     169        #endif
     170        return norm;
     171}
     172/*}}}*/
     173/*FUNCTION Matrix::GetSize{{{1*/
     174void Matrix::GetSize(int* pM,int* pN){
     175       
     176        #ifdef _HAVE_PETSC_
     177            MatGetSize(this->matrix,pM,pN);
     178        #else
     179                _error_("not implemented yet!");
     180        #endif
     181}
     182/*}}}*/
     183/*FUNCTION Matrix::GetLocalSize{{{1*/
     184void Matrix::GetLocalSize(int* pM,int* pN){
     185       
     186        #ifdef _HAVE_PETSC_
     187            MatGetLocalSize(this->matrix,pM,pN);
     188        #else
     189                _error_("not implemented yet!");
     190        #endif
     191}
     192/*}}}*/
     193/*FUNCTION Matrix::MatMult{{{1*/
     194void Matrix::MatMult(Vector* X,Vector* AX){
     195
     196        #ifdef _HAVE_PETSC_
     197            MatMultPatch(this->matrix,X->vector,AX->vector);
     198        #else
     199                _error_("not implemented yet!");
     200        #endif
     201}
     202/*}}}*/

  • issm/trunk-jpl/src/c/objects/Numerics/Matrix.h

    r11670 r11679  
    2525#include "mex.h"
    2626#endif
     27class Vector;
    2728
    2829/*}}}*/
    … …  
    5556                mxArray* ToMatlabMatrix(void);
    5657                #endif
     58                void Assemble(void);
     59                double Norm(int norm_type);
     60                void GetSize(int* pM,int* pN);
     61                void GetLocalSize(int* pM,int* pN);
     62                void MatMult(Vector* X,Vector* AX);
    5763                /*}}}*/
    5864

  • issm/trunk-jpl/src/c/objects/Numerics/Vector.cpp

    r11670 r11679  
    2424Vector::Vector(){
    2525
     26        this->M=0;
    2627        #ifdef _HAVE_PETSC_
    2728        this->vector=NULL;
    2829        #else
    2930        this->vector=NULL;
    30         this->M=0;
    3131        #endif
    3232        #ifdef _HAVE_ADOLC_
    … …  
    9898        _error_("not implemented yet!");
    9999        #endif
     100        return dataref;
    100101
    101102}
    102103/*}}}*/
    103104#endif
     105/*FUNCTION Vector::Assemble{{{1*/
     106void Vector::Assemble(void){
     107               
     108        #ifdef _HAVE_PETSC_
     109                VecAssemblyBegin(this->vector);
     110                VecAssemblyEnd(this->vector);
     111        #else
     112                /*do nothing*/
     113        #endif
     114
     115}
     116/*}}}*/
     117/*FUNCTION Vector::SetValues{{{1*/
     118void Vector::SetValues(int ssize, int* list, double* values, int mode){
     119               
     120               
     121        #ifdef _HAVE_PETSC_
     122        VecSetValues(this->vector,ssize,list,values,(InsertMode)mode);
     123        #else
     124                _error_("not implemented yet!");
     125        #endif
     126
     127}
     128/*}}}*/
     129/*FUNCTION Vector::GetSize{{{1*/
     130void Vector::GetSize(int* pM){
     131               
     132        #ifdef _HAVE_PETSC_
     133                VecGetSize(this->vector,pM);
     134        #else
     135                _error_("not implemented yet!");
     136        #endif
     137
     138}
     139/*}}}*/
     140/*FUNCTION Vector::GetLocalSize{{{1*/
     141void Vector::GetLocalSize(int* pM){
     142               
     143        #ifdef _HAVE_PETSC_
     144                VecGetLocalSize(this->vector,pM);
     145        #else
     146                _error_("not implemented yet!");
     147        #endif
     148
     149}
     150/*}}}*/
     151/*FUNCTION Vector::Duplicate{{{1*/
     152Vector* Vector::Duplicate(void){
     153       
     154        Vector* output=NULL;
     155        output=new Vector();
     156        #ifdef _HAVE_PETSC_
     157                Vec vec_output=NULL;
     158                VecDuplicate(this->vector,&vec_output);
     159                output->vector=vec_output;
     160                VecGetSize(output->vector,&output->M);
     161        #else
     162                _error_("not implemented yet!");
     163        #endif
     164        return output;
     165
     166}
     167/*}}}*/
     168/*FUNCTION Vector::Set{{{1*/
     169void Vector::Set(double value){
     170       
     171        #ifdef _HAVE_PETSC_
     172                this->Set(value);
     173        #else
     174                _error_("not implemented yet!");
     175        #endif
     176
     177}
     178/*}}}*/
     179/*FUNCTION Vector::AXPY{{{1*/
     180void Vector::AXPY(Vector* X, double a){
     181       
     182        #ifdef _HAVE_PETSC_
     183                VecAXPY(this->vector,a,X->vector);
     184        #else
     185                _error_("not implemented yet!");
     186        #endif
     187
     188}
     189/*}}}*/
     190/*FUNCTION Vector::AYPX{{{1*/
     191void Vector::AYPX(Vector* X, double a){
     192       
     193        #ifdef _HAVE_PETSC_
     194                VecAYPX(this->vector,a,X->vector);
     195        #else
     196                _error_("not implemented yet!");
     197        #endif
     198
     199}
     200/*}}}*/
     201/*FUNCTION Vector::ToMPISerial{{{1*/
     202double* Vector::ToMPISerial(void){
     203
     204        double* vec_serial=NULL;
     205
     206        #ifdef _HAVE_PETSC_
     207                VecToMPISerial(&vec_serial, this->vector);
     208        #else
     209                _error_("not implemented yet!");
     210        #endif
     211
     212        return vec_serial;
     213
     214}
     215/*}}}*/
     216/*FUNCTION Vector::Copy{{{1*/
     217void Vector::Copy(Vector* to){
     218
     219        #ifdef _HAVE_PETSC_
     220                VecCopy(this->vector,to->vector);
     221        #else
     222                _error_("not implemented yet!");
     223        #endif
     224
     225}
     226/*}}}*/
     227/*FUNCTION Vector::Norm{{{1*/
     228double Vector::Norm(int norm_type){
     229       
     230        double norm=0;
     231        #ifdef _HAVE_PETSC_
     232            VecNorm(this->vector,(NormType)norm_type,&norm);
     233        #else
     234                _error_("not implemented yet!");
     235        #endif
     236        return norm;
     237}
     238/*}}}*/
     239/*FUNCTION Vector::Scale{{{1*/
     240void Vector::Scale(double scale_factor){
     241       
     242        #ifdef _HAVE_PETSC_
     243            VecScale(this->vector,scale_factor);
     244        #else
     245                _error_("not implemented yet!");
     246        #endif
     247}
     248/*}}}*/

  • issm/trunk-jpl/src/c/objects/Numerics/Vector.h

    r11670 r11679  
    5454                mxArray* ToMatlabVector(void);
    5555                #endif
     56                void Assemble(void);
     57                void SetValues(int ssize, int* list, double* values, int mode);
     58                void GetSize(int* pM);
     59                void GetLocalSize(int* pM);
     60                Vector* Duplicate(void);
     61                void Set(double value);
     62                void AXPY(Vector* X, double a);
     63                void AYPX(Vector* X, double a);
     64                double* ToMPISerial(void);
     65                void Copy(Vector* to);
     66                double Norm(int norm_type);
     67                void Scale(double scale_factor);
    5668                /*}}}*/
    5769};

  • issm/trunk-jpl/src/c/shared/Alloc/alloc.cpp

    r9320 r11679  
    2929#include "../Exceptions/exceptions.h"
    3030#include "../../include/include.h"
     31#include "../../objects/objects.h"
    3132
    3233void* xmalloc(int size){
    … …  
    8081}
    8182
     83void xdelete( Matrix** pv){
     84
     85        if (pv && *pv) {
     86
     87                delete *pv;
     88                *pv=NULL;
     89        }
     90}
     91
     92void xdelete( Vector** pv){
     93
     94        if (pv && *pv) {
     95
     96                delete *pv;
     97                *pv=NULL;
     98        }
     99}
     100
     101
    82102void* xrealloc ( void* pv, int size){
    83103       

  • issm/trunk-jpl/src/c/shared/Alloc/alloc.h

    r5768 r11679  
    55#ifndef _ALLOC_H_
    66#define _ALLOC_H_
    7 
     7class Matrix;
     8class Vector;
    89void* xmalloc(int size);
    910void* xcalloc(int n,int size);
    1011void  xfree(void** pvptr);
    1112void* xrealloc ( void* pv, int size);
     13void xdelete(Matrix** pvptr);
     14void xdelete(Vector** pvptr);
    1215
    1316#endif

  • issm/trunk-jpl/src/c/solutions/ResetBoundaryConditions.cpp

    r9761 r11679  
    1111       
    1212        /*variables: */
    13         Vec    yg    = NULL;
     13        Vector*    yg    = NULL;
    1414        Nodes *nodes = NULL;
    1515        int    i;
    … …  
    3030
    3131        /*Free ressources:*/
    32         VecFree(&yg);
     32        delete yg;
    3333}

  • issm/trunk-jpl/src/c/solutions/convergence.cpp

    r11285 r11679  
    88#include "../EnumDefinitions/EnumDefinitions.h"
    99
    10 void convergence(bool* pconverged, Mat Kff,Vec pf,Vec uf,Vec old_uf,Parameters* parameters){
     10void convergence(bool* pconverged, Matrix* Kff,Vector* pf,Vector* uf,Vector* old_uf,Parameters* parameters){
    1111
    1212        /*output*/
    … …  
    1414
    1515        /*intermediary*/
    16         Vec KU=NULL;
    17         Vec KUF=NULL;
    18         Vec KUold=NULL;
    19         Vec KUoldF=NULL;
    20         Vec duf=NULL;
     16        Vector* KU=NULL;
     17        Vector* KUF=NULL;
     18        Vector* KUold=NULL;
     19        Vector* KUoldF=NULL;
     20        Vector* duf=NULL;
    2121        double ndu,nduinf,nu;
    2222        double nKUF;
    … …  
    4848
    4949                //compute KUF = KU - F = K*U - F
    50                 VecDuplicate(uf,&KU); MatMultPatch(Kff,uf,KU);
    51                 VecDuplicate(KU,&KUF);VecCopy(KU,KUF); VecAYPX(KUF,-1.0,pf);
     50                KU=uf->Duplicate(); Kff->MatMult(uf,KU);
     51                KUF=KU->Duplicate(); KU->Copy(KUF); KUF->AYPX(pf,-1.0);
    5252
    5353                //compute norm(KUF), norm(F) and residue
    54                 VecNorm(KUF,NORM_2,&nKUF);
    55                 VecNorm(pf,NORM_2,&nF);
     54                nKUF=KUF->Norm(NORM_2);
     55                nF=pf->Norm(NORM_2);
    5656                solver_residue=nKUF/nF;
    5757                _printf_(true,"\n%s%g\n","   solver residue: norm(KU-F)/norm(F)=",solver_residue);
    5858
    5959                //clean up
    60                 VecFree(&KU);
    61                 VecFree(&KUF);
     60                delete KU;
     61                delete KUF;
    6262        }
    6363
    … …  
    6565
    6666        //compute K[n]U[n-1]F = K[n]U[n-1] - F
    67         VecDuplicate(uf,&KUold); MatMultPatch(Kff,old_uf,KUold);
    68         VecDuplicate(KUold,&KUoldF);VecCopy(KUold,KUoldF); VecAYPX(KUoldF,-1.0,pf);
    69         VecNorm(KUoldF,NORM_2,&nKUoldF);
    70         VecNorm(pf,NORM_2,&nF);
     67        KUold=uf->Duplicate(); Kff->MatMult(old_uf,KUold);
     68        KUoldF=KUold->Duplicate();KUold->Copy(KUoldF); KUoldF->AYPX(pf,-1.0);
     69        nKUoldF=KUoldF->Norm(NORM_2);
     70        nF=pf->Norm(NORM_2);
    7171        res=nKUoldF/nF;
    7272        if (isnan(res)){
    … …  
    7676
    7777        //clean up
    78         VecFree(&KUold);
    79         VecFree(&KUoldF);
     78        delete KUold;
     79        delete KUoldF;
    8080
    8181        //print
    … …  
    9393
    9494                //compute norm(du)/norm(u)
    95                 VecDuplicate(old_uf,&duf);VecCopy(old_uf,duf); VecAYPX(duf,-1.0,uf);
    96                 VecNorm(duf,NORM_2,&ndu); VecNorm(old_uf,NORM_2,&nu);
     95                duf=old_uf->Duplicate(); old_uf->Copy(duf); duf->AYPX(uf,-1.0);
     96                ndu=duf->Norm(NORM_2); nu=old_uf->Norm(NORM_2);
    9797
    9898                if (isnan(ndu) || isnan(nu)) _error_("convergence criterion is NaN!");
    9999
    100100                //clean up
    101                 VecFree(&duf);
     101                delete duf;
    102102
    103103                //print
    … …  
    119119
    120120                //compute max(du)
    121                 VecDuplicate(old_uf,&duf);VecCopy(old_uf,duf); VecAYPX(duf,-1.0,uf);
    122                 VecNorm(duf,NORM_2,&ndu); VecNorm(duf,NORM_INFINITY,&nduinf);
     121                duf=old_uf->Duplicate(); old_uf->Copy(duf); duf->AYPX(uf,-1.0);
     122                ndu=duf->Norm(NORM_2); nduinf=duf->Norm(NORM_INFINITY);
    123123                if (isnan(ndu) || isnan(nu)) _error_("convergence criterion is NaN!");
    124124
    125125                //clean up
    126                 VecFree(&duf);
     126                delete duf;
    127127
    128128                //print

  • issm/trunk-jpl/src/c/solutions/solutions.h

    r11306 r11679  
    3535
    3636//convergence:
    37 void convergence(bool* pconverged, Mat K_ff,Vec p_f,Vec u_f,Vec u_f_old,Parameters* parameters);
     37void convergence(bool* pconverged, Matrix* K_ff,Vector* p_f,Vector* u_f,Vector* u_f_old,Parameters* parameters);
    3838bool controlconvergence(double J,double tol_cm);
    3939bool steadystateconvergence(FemModel* femmodel);

  • issm/trunk-jpl/src/c/solvers/solver_adjoint_linear.cpp

    r9271 r11679  
    1313
    1414        /*intermediary: */
    15         Mat  Kff = NULL, Kfs = NULL;
    16         Vec  ug  = NULL, uf  = NULL;
    17         Vec  pf  = NULL;
    18         Vec  df  = NULL;
    19         Vec  ys  = NULL;
     15        Matrix*  Kff = NULL;
     16        Matrix*  Kfs = NULL;
     17        Vector*  ug  = NULL;
     18        Vector*  uf  = NULL;
     19        Vector*  pf  = NULL;
     20        Vector*  df  = NULL;
     21        Vector*  ys  = NULL;
    2022        int  configuration_type;
    2123
    … …  
    2628        SystemMatricesx(&Kff, &Kfs, &pf, &df, NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters);
    2729        CreateNodalConstraintsx(&ys,femmodel->nodes,configuration_type);
    28         Reduceloadx(pf, Kfs, ys,true); MatFree(&Kfs); //true means spc = 0
    29         Solverx(&uf, Kff, pf, NULL, df, femmodel->parameters); MatFree(&Kff); VecFree(&pf); VecFree(&df);
    30         Mergesolutionfromftogx(&ug, uf,ys,femmodel->nodes,femmodel->parameters,true); VecFree(&uf);VecFree(&ys); //true means spc0
     30        Reduceloadx(pf, Kfs, ys,true); xdelete(&Kfs); //true means spc = 0
     31        Solverx(&uf, Kff, pf, NULL, df, femmodel->parameters); xdelete(&Kff); xdelete(&pf); xdelete(&df);
     32        Mergesolutionfromftogx(&ug, uf,ys,femmodel->nodes,femmodel->parameters,true); xdelete(&uf);xdelete(&ys); //true means spc0
    3133        InputUpdateFromSolutionx(femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,ug);
    32         VecFree(&ug); VecFree(&uf);
     34        xdelete(&ug); xdelete(&uf);
    3335}

  • issm/trunk-jpl/src/c/solvers/solver_linear.cpp

    r9225 r11679  
    1111
    1212        /*intermediary: */
    13         Mat  Kff = NULL, Kfs   = NULL;
    14         Vec  ug  = NULL;
    15         Vec  uf  = NULL;
    16         Vec  pf  = NULL;
    17         Vec  df  = NULL;
    18         Vec  ys  = NULL;
     13        Matrix*  Kff = NULL;
     14        Matrix*  Kfs = NULL;
     15        Vector*  ug  = NULL;
     16        Vector*  uf  = NULL;
     17        Vector*  pf  = NULL;
     18        Vector*  df  = NULL;
     19        Vector*  ys  = NULL;
    1920        int  configuration_type;
    2021
    … …  
    2526        SystemMatricesx(&Kff, &Kfs, &pf, &df, NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters);
    2627        CreateNodalConstraintsx(&ys,femmodel->nodes,configuration_type);
    27         Reduceloadx(pf, Kfs, ys); MatFree(&Kfs);
    28         Solverx(&uf, Kff, pf, NULL, df, femmodel->parameters); MatFree(&Kff); VecFree(&pf); VecFree(&df);
    29         Mergesolutionfromftogx(&ug, uf,ys,femmodel->nodes,femmodel->parameters);VecFree(&uf);VecFree(&ys);
     28        Reduceloadx(pf, Kfs, ys); xdelete(&Kfs);
     29        Solverx(&uf, Kff, pf, NULL, df, femmodel->parameters);
     30        xdelete(&Kff); xdelete(&pf); xdelete(&df);
     31        Mergesolutionfromftogx(&ug, uf,ys,femmodel->nodes,femmodel->parameters);xdelete(&uf); xdelete(&ys);
    3032        InputUpdateFromSolutionx( femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,ug);
    31         VecFree(&ug); VecFree(&uf);
     33        xdelete(&ug);  xdelete(&uf);
    3234}

  • issm/trunk-jpl/src/c/solvers/solver_newton.cpp

    r11338 r11679  
    1818        int    count;
    1919        double kmax;
    20         Mat Kff = NULL, Kfs    = NULL, Jff = NULL;
    21         Vec ug  = NULL, old_ug = NULL;
    22         Vec uf  = NULL, old_uf = NULL, duf = NULL;
    23         Vec pf  = NULL, pJf    = NULL;
    24         Vec df  = NULL;
    25         Vec ys  = NULL;
     20        Matrix* Kff = NULL;
     21        Matrix* Kfs    = NULL;
     22        Matrix* Jff = NULL;
     23        Vector* ug  = NULL;
     24        Vector* old_ug = NULL;
     25        Vector* uf  = NULL;
     26        Vector* old_uf = NULL;
     27        Vector* duf = NULL;
     28        Vector* pf  = NULL;
     29        Vector* pJf    = NULL;
     30        Vector* df  = NULL;
     31        Vector* ys  = NULL;
    2632
    2733        /*parameters:*/
    … …  
    4753        for(;;){
    4854
    49                 VecFree(&old_ug);old_ug=ug;
    50                 VecFree(&old_uf);old_uf=uf;
     55                xdelete(&old_ug);old_ug=ug;
     56                xdelete(&old_uf);old_uf=uf;
    5157
    5258                /*Solver forward model*/
    5359                SystemMatricesx(&Kff,&Kfs,&pf,&df,NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters);
    5460                CreateNodalConstraintsx(&ys,femmodel->nodes,configuration_type);
    55                 Reduceloadx(pf,Kfs,ys);MatFree(&Kfs);
    56                 Solverx(&uf,Kff,pf,old_uf,df,femmodel->parameters);VecFree(&df);
    57                 Mergesolutionfromftogx(&ug,uf,ys,femmodel->nodes,femmodel->parameters);VecFree(&ys);
    58                 InputUpdateFromSolutionx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,ug);VecFree(&ug);
     61                Reduceloadx(pf,Kfs,ys);xdelete(&Kfs);
     62                Solverx(&uf,Kff,pf,old_uf,df,femmodel->parameters);xdelete(&df);
     63                Mergesolutionfromftogx(&ug,uf,ys,femmodel->nodes,femmodel->parameters);xdelete(&ys);
     64                InputUpdateFromSolutionx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,ug);xdelete(&ug);
    5965
    6066                /*Check convergence*/
    6167                convergence(&converged,Kff,pf,uf,old_uf,femmodel->parameters);
    62                 MatFree(&Kff);VecFree(&pf);
     68                xdelete(&Kff); xdelete(&pf);
    6369                if(converged==true) break;
    6470                if(count>=max_nonlinear_iterations){
    … …  
    7076                SystemMatricesx(&Kff,&Kfs,&pf,NULL,&kmax,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters);
    7177                CreateNodalConstraintsx(&ys,femmodel->nodes,configuration_type);
    72                 Reduceloadx(pf,Kfs,ys);   MatFree(&Kfs);
     78                Reduceloadx(pf,Kfs,ys);   xdelete(&Kfs);
    7379
    74                 VecDuplicate(pf,&pJf);
    75                 MatMultPatch(Kff,uf,pJf); MatFree(&Kff);
    76                 VecScale(pJf,-1.);
    77                 VecAXPY(pJf,+1.,pf);      VecFree(&pf);
     80                pJf=pf->Duplicate(); Kff->MatMult(uf,pJf); xdelete(&Kff);
     81                pJf->Scale(-1.0); pJf->AXPY(pf,+1.0);     xdelete(&pf);
    7882
    7983                CreateJacobianMatrixx(&Jff,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,kmax);
    80                 Solverx(&duf,Jff,pJf,NULL,NULL,femmodel->parameters); MatFree(&Jff);VecFree(&pJf);
    81                 VecAXPY(uf,1.,duf);      VecFree(&duf);
    82                 Mergesolutionfromftogx(&ug,uf,ys,femmodel->nodes,femmodel->parameters);VecFree(&ys);
     84                Solverx(&duf,Jff,pJf,NULL,NULL,femmodel->parameters); xdelete(&Jff); xdelete(&pJf);
     85                uf->AXPY(duf, 1.0); xdelete(&duf);
     86                Mergesolutionfromftogx(&ug,uf,ys,femmodel->nodes,femmodel->parameters);xdelete(&ys);
    8387                InputUpdateFromSolutionx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,ug);
    8488
    … …  
    8993
    9094        /*clean-up*/
    91         VecFree(&uf);
    92         VecFree(&ug);
    93         VecFree(&old_ug);
    94         VecFree(&old_uf);
     95        xdelete(&uf);
     96        xdelete(&ug);
     97        xdelete(&old_ug);
     98        xdelete(&old_uf);
    9599}

  • issm/trunk-jpl/src/c/solvers/solver_nonlinear.cpp

    r11347 r11679  
    1414
    1515        /*intermediary: */
    16         Mat Kff = NULL, Kfs   = NULL;
    17         Vec ug  = NULL, old_ug = NULL;
    18         Vec uf  = NULL, old_uf = NULL;
    19         Vec pf  = NULL;
    20         Vec df  = NULL;
    21         Vec ys  = NULL;
     16        Matrix* Kff = NULL;
     17        Matrix* Kfs = NULL;
     18        Vector* ug  = NULL;
     19        Vector* old_ug = NULL;
     20        Vector* uf  = NULL;
     21        Vector* old_uf = NULL;
     22        Vector* pf  = NULL;
     23        Vector* df  = NULL;
     24        Vector* ys  = NULL;
    2225       
    2326        Loads* loads=NULL;
    … …  
    5659
    5760                //save pointer to old velocity
    58                 VecFree(&old_ug);old_ug=ug;
    59                 VecFree(&old_uf);old_uf=uf;
     61                xdelete(&old_ug);old_ug=ug;
     62                xdelete(&old_uf);old_uf=uf;
    6063
    6164                SystemMatricesx(&Kff, &Kfs, &pf, &df, NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,loads,femmodel->materials,femmodel->parameters);
    6265                CreateNodalConstraintsx(&ys,femmodel->nodes,configuration_type);
    63                 Reduceloadx(pf, Kfs, ys); MatFree(&Kfs);
     66                Reduceloadx(pf, Kfs, ys); xdelete(&Kfs);
    6467                Solverx(&uf, Kff, pf, old_uf, df, femmodel->parameters);
    65                 Mergesolutionfromftogx(&ug, uf,ys,femmodel->nodes,femmodel->parameters);VecFree(&ys);
     68                Mergesolutionfromftogx(&ug, uf,ys,femmodel->nodes,femmodel->parameters);xdelete(&ys);
    6669
    67                 convergence(&converged,Kff,pf,uf,old_uf,femmodel->parameters); MatFree(&Kff);VecFree(&pf); VecFree(&df);
     70                convergence(&converged,Kff,pf,uf,old_uf,femmodel->parameters); xdelete(&Kff); xdelete(&pf); xdelete(&df);
    6871                InputUpdateFromConstantx( femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,converged,ConvergedEnum);
    6972                InputUpdateFromSolutionx( femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,ug);
    … …  
    9699        /*clean-up*/
    97100        if(conserve_loads) delete loads;
    98         VecFree(&uf);
    99         VecFree(&ug);
    100         VecFree(&old_ug);
    101         VecFree(&old_uf);
     101        xdelete(&uf);
     102        xdelete(&ug);
     103        xdelete(&old_ug);
     104        xdelete(&old_uf);
    102105}

  • issm/trunk-jpl/src/c/solvers/solver_stokescoupling_nonlinear.cpp

    r11284 r11679  
    1414
    1515        /*intermediary: */
    16         Mat  Kff_horiz = NULL, Kfs_horiz   = NULL;
    17         Vec  ug_horiz  = NULL, uf_horiz  = NULL, old_uf_horiz = NULL;
    18         Vec  pf_horiz  = NULL;
    19         Vec  df_horiz  = NULL;
    20         Mat  Kff_vert  = NULL, Kfs_vert    = NULL;
    21         Vec  ug_vert   = NULL, uf_vert   = NULL;
    22         Vec  pf_vert   = NULL;
    23         Vec  df_vert   = NULL;
    24         Vec  ys   = NULL;
     16        Matrix*  Kff_horiz = NULL;
     17        Matrix* Kfs_horiz   = NULL;
     18        Vector*  ug_horiz  = NULL;
     19        Vector*  uf_horiz  = NULL;
     20        Vector*  old_uf_horiz = NULL;
     21        Vector*  pf_horiz  = NULL;
     22        Vector*  df_horiz  = NULL;
     23        Matrix*  Kff_vert  = NULL;
     24        Matrix*  Kfs_vert  = NULL;
     25        Vector*  ug_vert   = NULL;
     26        Vector*  uf_vert   = NULL;
     27        Vector*  pf_vert   = NULL;
     28        Vector*  df_vert   = NULL;
     29        Vector*  ys   = NULL;
    2530        bool converged;
    2631        int  constraints_converged;
    … …  
    5459                //Update once again the solution to make sure that vx and vxold are similar (for next step in transient or steadystate)
    5560                InputUpdateFromSolutionx(femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,ug_horiz);
    56                 VecFree(&ug_horiz);
     61                xdelete(&ug_horiz);
    5762
    5863                //save pointer to old velocity
    59                 VecFree(&old_uf_horiz);old_uf_horiz=uf_horiz;
     64                xdelete(&old_uf_horiz); old_uf_horiz=uf_horiz;
    6065
    6166                /*solve: */
    6267                SystemMatricesx(&Kff_horiz, &Kfs_horiz, &pf_horiz, &df_horiz, NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters);
    6368                CreateNodalConstraintsx(&ys,femmodel->nodes,configuration_type);
    64                 Reduceloadx(pf_horiz, Kfs_horiz, ys); MatFree(&Kfs_horiz);
     69                Reduceloadx(pf_horiz, Kfs_horiz, ys); xdelete(&Kfs_horiz);
    6570                Solverx(&uf_horiz, Kff_horiz, pf_horiz, old_uf_horiz, df_horiz,femmodel->parameters);
    66                 Mergesolutionfromftogx(&ug_horiz, uf_horiz,ys,femmodel->nodes,femmodel->parameters); VecFree(&ys);
     71                Mergesolutionfromftogx(&ug_horiz, uf_horiz,ys,femmodel->nodes,femmodel->parameters); xdelete(&ys);
    6772                InputUpdateFromSolutionx( femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,ug_horiz);
    6873
    69                 convergence(&converged,Kff_horiz,pf_horiz,uf_horiz,old_uf_horiz,femmodel->parameters); MatFree(&Kff_horiz);VecFree(&pf_horiz); VecFree(&df_horiz);
     74                convergence(&converged,Kff_horiz,pf_horiz,uf_horiz,old_uf_horiz,femmodel->parameters); xdelete(&Kff_horiz); xdelete(&pf_horiz); xdelete(&df_horiz);
    7075
    7176                /*Second compute vertical velocity: */
    … …  
    7681                SystemMatricesx(&Kff_vert, &Kfs_vert, &pf_vert,  &df_vert,NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters);
    7782                CreateNodalConstraintsx(&ys,femmodel->nodes,configuration_type);
    78                 Reduceloadx(pf_vert, Kfs_vert, ys); MatFree(&Kfs_vert);
    79                 Solverx(&uf_vert, Kff_vert, pf_vert, NULL, df_vert,femmodel->parameters); MatFree(&Kff_vert); VecFree(&pf_vert); VecFree(&df_vert);
    80                 Mergesolutionfromftogx(&ug_vert, uf_vert,ys,femmodel->nodes,femmodel->parameters);VecFree(&uf_vert); VecFree(&ys);
     83                Reduceloadx(pf_vert, Kfs_vert, ys); xdelete(&Kfs_vert);
     84                Solverx(&uf_vert, Kff_vert, pf_vert, NULL, df_vert,femmodel->parameters); xdelete(&Kff_vert); xdelete(&pf_vert); xdelete(&df_vert);
     85                Mergesolutionfromftogx(&ug_vert, uf_vert,ys,femmodel->nodes,femmodel->parameters);xdelete(&uf_vert); xdelete(&ys);
    8186                InputUpdateFromSolutionx( femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,ug_vert);
    82                 VecFree(&ug_vert); VecFree(&uf_vert);
     87                xdelete(&ug_vert); xdelete(&uf_vert);
    8388
    8489                /*Increase count: */
    … …  
    9297
    9398        /*clean-up*/
    94         VecFree(&old_uf_horiz);
    95         VecFree(&uf_horiz);
    96         VecFree(&ug_horiz);
    97         VecFree(&ys);
     99        xdelete(&old_uf_horiz);
     100        xdelete(&uf_horiz);
     101        xdelete(&ug_horiz);
     102        xdelete(&ys);
    98103}

  • issm/trunk-jpl/src/c/solvers/solver_thermal_nonlinear.cpp

    r11197 r11679  
    1212
    1313        /*solution : */
    14         Vec tg=NULL;
    15         Vec tf=NULL;
    16         Vec tf_old=NULL;
    17         Vec ys=NULL;
     14        Vector* tg=NULL;
     15        Vector* tf=NULL;
     16        Vector* tf_old=NULL;
     17        Vector* ys=NULL;
    1818        double melting_offset;
    1919
    2020        /*intermediary: */
    21         Mat Kff=NULL;
    22         Mat Kfs=NULL;
    23         Vec pf=NULL;
    24         Vec df=NULL;
     21        Matrix* Kff=NULL;
     22        Matrix* Kfs=NULL;
     23        Vector* pf=NULL;
     24        Vector* df=NULL;
    2525
    2626        bool converged;
    … …  
    5656                SystemMatricesx(&Kff, &Kfs, &pf,&df, &melting_offset,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters);
    5757                CreateNodalConstraintsx(&ys,femmodel->nodes,configuration_type);
    58                 Reduceloadx(pf, Kfs, ys); MatFree(&Kfs); VecFree(&tf);
     58                Reduceloadx(pf, Kfs, ys); xdelete(&Kfs); xdelete(&tf);
    5959                Solverx(&tf, Kff, pf,tf_old, df, femmodel->parameters);
    60                 VecFree(&tf_old); VecDuplicatePatch(&tf_old,tf);
    61                 MatFree(&Kff);VecFree(&pf);VecFree(&tg); VecFree(&df);
    62                 Mergesolutionfromftogx(&tg, tf,ys,femmodel->nodes,femmodel->parameters); VecFree(&ys);
     60                xdelete(&tf_old); tf_old=tf->Duplicate();
     61                xdelete(&Kff);xdelete(&pf);xdelete(&tg); xdelete(&df);
     62                Mergesolutionfromftogx(&tg, tf,ys,femmodel->nodes,femmodel->parameters); xdelete(&ys);
    6363                InputUpdateFromSolutionx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,tg);
    6464
    … …  
    8484
    8585        /*Free ressources: */
    86         VecFree(&tg);
    87         VecFree(&tf);
    88         VecFree(&tf_old);
    89         VecFree(&ys);
     86        xdelete(&tg);
     87        xdelete(&tf);
     88        xdelete(&tf_old);
     89        xdelete(&ys);
    9090}

  • issm/trunk-jpl/src/mex/CreateJacobianMatrix/CreateJacobianMatrix.cpp

    r11335 r11679  
    1717       
    1818        /* output datasets: */
    19         Mat    Jff  = NULL;
     19        Matrix*    Jff  = NULL;
    2020
    2121        /*Boot module: */
    … …  
    5353        delete materials;
    5454        delete parameters;
    55         MatFree(&Jff);
     55        delete Jff;
    5656
    5757        /*end module: */

  • issm/trunk-jpl/src/mex/CreateNodalConstraints/CreateNodalConstraints.cpp

    r8910 r11679  
    1212
    1313        /* output datasets: */
    14         Vec ys=NULL;
     14        Vector* ys=NULL;
    1515
    1616        /*Boot module: */
    … …  
    3232        /*Free ressources: */
    3333        delete nodes;
    34         VecFree(&ys);
     34        delete ys;
    3535
    3636        /*end module: */

  • issm/trunk-jpl/src/mex/GetSolutionFromInputs/GetSolutionFromInputs.cpp

    r8910 r11679  
    1414        Materials* materials=NULL;
    1515        Parameters* parameters=NULL;
    16         Vec      ug=NULL;
     16        Vector*      ug=NULL;
    1717
    1818        /* output datasets: elements and loads*/

  • issm/trunk-jpl/src/mex/InputUpdateFromSolution/InputUpdateFromSolution.cpp

    r8910 r11679  
    1414        Materials* materials=NULL;
    1515        Parameters* parameters=NULL;
    16         Vec      solution=NULL;
     16        Vector*      solution=NULL;
    1717
    1818        /*Boot module: */
    … …  
    5050        delete materials;
    5151        delete parameters;
    52         VecFree(&solution);
     52        delete solution;
    5353
    5454        /*end module: */

  • issm/trunk-jpl/src/mex/Mergesolutionfromftog/Mergesolutionfromftog.cpp

    r8910 r11679  
    99        /*input datasets: */
    1010        bool        flag_ys0;
    11         Vec         uf         = NULL;
    12         Vec         ys         = NULL;
     11        Vector*         uf         = NULL;
     12        Vector*         ys         = NULL;
    1313        Nodes*      nodes   = NULL;
    1414        Parameters* parameters   = NULL;
    1515
    1616        /* output datasets: */
    17         Vec ug=NULL;
     17        Vector* ug=NULL;
    1818
    1919        /*Boot module: */
    … …  
    4545
    4646        /*Free ressources: */
    47         VecFree(&uf);
    48         VecFree(&ug);
    49         VecFree(&ys);
     47        delete uf;
     48        delete ug;
     49        delete ys;
    5050        delete nodes;
    5151        delete parameters;

  • issm/trunk-jpl/src/mex/Reduceload/Reduceload.cpp

    r8910 r11679  
    88
    99        /*input datasets: */
    10         Vec         pf         = NULL;
    11         Mat         Kfs        = NULL;
    12         Vec         ys         = NULL;
     10        Vector*         pf         = NULL;
     11        Matrix*         Kfs        = NULL;
     12        Vector*         ys         = NULL;
    1313        bool        flag_ys0=false;
    1414
    … …  
    4040
    4141        /*Free ressources: */
    42         VecFree(&pf);
    43         MatFree(&Kfs);
    44         VecFree(&ys);
     42        delete pf;
     43        delete Kfs;
     44        delete ys;
    4545
    4646        MODULEEND();

  • issm/trunk-jpl/src/mex/Reducevectorgtof/Reducevectorgtof.cpp

    r8910 r11679  
    88
    99        /*input datasets: */
    10         Vec ug=NULL;
     10        Vector* ug=NULL;
    1111        Nodes* nodes=NULL;
    1212        Parameters* parameters=NULL;
    1313
    1414        /* output datasets: */
    15         Vec uf=NULL;
     15        Vector* uf=NULL;
    1616
    1717        /*Boot module: */
    … …  
    3535        delete nodes;
    3636        delete parameters;
    37         VecFree(&ug);
    38         VecFree(&uf);
     37        delete ug;
     38        delete uf;
    3939
    4040        /*end module: */

  • issm/trunk-jpl/src/mex/Reducevectorgtos/Reducevectorgtos.cpp

    r8910 r11679  
    88
    99        /*input datasets: */
    10         Vec yg=NULL;
     10        Vector* yg=NULL;
    1111        Nodes* nodes=NULL;
    1212        Parameters* parameters=NULL;
    1313
    1414        /* output datasets: */
    15         Vec ys=NULL;
     15        Vector* ys=NULL;
    1616
    1717        /*Boot module: */
    … …  
    3535        delete nodes;
    3636        delete parameters;
    37         VecFree(&yg);
    38         VecFree(&ys);
     37        delete yg;
     38        delete ys;
    3939
    4040        /*end module: */

  • issm/trunk-jpl/src/mex/Solver/Solver.cpp

    r11252 r11679  
    88
    99        /*input datasets: */
    10         Mat         Kff           = NULL;
    11         Vec         pf            = NULL;
    12         Vec         uf0           = NULL;
    13         Vec         uf            = NULL;
    14         Vec         df            = NULL;
     10        Matrix*         Kff           = NULL;
     11        Vector*         pf            = NULL;
     12        Vector*         uf0           = NULL;
     13        Vector*         uf            = NULL;
     14        Vector*         df            = NULL;
    1515        Parameters *parameters    = NULL;
    1616        int         analysis_type;
    … …  
    6666
    6767        /*Free ressources: */
    68         MatFree(&Kff);
    69         VecFree(&pf);
    70         VecFree(&uf0);
    71         VecFree(&uf);
    72         VecFree(&df);
     68        delete Kff;
     69        delete pf;
     70        delete uf0;
     71        delete uf;
     72        delete df;
    7373        delete parameters;
    7474

  • issm/trunk-jpl/src/mex/SystemMatrices/SystemMatrices.cpp

    r8910 r11679  
    1717       
    1818        /* output datasets: */
    19         Mat    Kff  = NULL;
    20         Mat    Kfs  = NULL;
    21         Vec    pf   = NULL;
    22         Vec    df   = NULL;
     19        Matrix*    Kff  = NULL;
     20        Matrix*    Kfs  = NULL;
     21        Vector*    pf   = NULL;
     22        Vector*    df   = NULL;
    2323
    2424        double kmax;
    … …  
    7272        delete materials;
    7373        delete parameters;
    74         MatFree(&Kff);
    75         MatFree(&Kfs);
    76         VecFree(&pf);
    77         VecFree(&df);
     74        delete Kff;
     75        delete Kfs;
     76        delete pf;
     77        delete df;
    7878
    7979        /*end module: */

  • issm/trunk-jpl/src/mex/UpdateDynamicConstraints/UpdateDynamicConstraints.cpp

    r9302 r11679  
    1111        Nodes       *nodes       = NULL;
    1212        Parameters  *parameters  = NULL;
    13         Vec          yg          = NULL;
     13        Vector*          yg          = NULL;
    1414
    1515        /*Boot module: */
    … …  
    3232
    3333        /*Free ressources: */
    34         VecFree(&yg);
     34        delete yg;
    3535        delete constraints;
    3636        delete nodes;
Note: See TracChangeset for help on using the changeset viewer.