source: issm/oecreview/Archive/11730-11748/ISSM-11733-11734.diff

Last change on this file was 11991, checked in by Eric.Larour, 13 years ago

oecreview from 11518 to present
File size: 20.1 KB

  • proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Numerics/ElementMatrix.cpp

     
    259259        /*In debugging mode, check consistency (no NaN, and values not too big)*/
    260260        this->CheckConsistency();
    261261
    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                                         }
    273                                 }
    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);
     262        if(this->dofsymmetrical){
     263                /*only use row dofs to add values into global matrices: */
    276264
    277                                 /*Free ressources:*/
    278                                 xfree((void**)&localvalues);
     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]);
     271                                }
    279272                        }
     273                        /*add local values into global  matrix, using the fglobaldoflist: */
     274                        Kff->SetValues(this->row_fsize,this->row_fglobaldoflist,this->row_fsize,this->row_fglobaldoflist,localvalues,ADD_VAL);
    280275
     276                        /*Free ressources:*/
     277                        xfree((void**)&localvalues);
     278                }
    281279
    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                                         }
    289                                 }
    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);
    292280
    293                                 /*Free ressources:*/
    294                                 xfree((void**)&localvalues);
     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]);
     287                                }
    295288                        }
     289                        /*add local values into global  matrix, using the fglobaldoflist: */
     290                        Kfs->SetValues(this->row_fsize,this->row_fglobaldoflist,this->row_ssize,this->row_sglobaldoflist,localvalues,ADD_VAL);
     291
     292                        /*Free ressources:*/
     293                        xfree((void**)&localvalues);
    296294                }
    297                 else{
    298                         _error_(" non dofsymmetrical matrix AddToGlobal routine not support yet!");
    299                 }
    300         #else
    301                 _error_("not supported yet!");
    302         #endif
     295        }
     296        else{
     297                _error_(" non dofsymmetrical matrix AddToGlobal routine not support yet!");
     298        }
    303299
    304300}
    305301/*}}}*/
    … …  
    315311        /*In debugging mode, check consistency (no NaN, and values not too big)*/
    316312        this->CheckConsistency();
    317313
    318         #ifdef _HAVE_PETSC_
    319                 if(this->dofsymmetrical){
    320                         /*only use row dofs to add values into global matrices: */
     314        if(this->dofsymmetrical){
     315                /*only use row dofs to add values into global matrices: */
    321316
    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                                         }
     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]);
    329323                                }
    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);
    335324                        }
     325                        /*add local values into global  matrix, using the fglobaldoflist: */
     326                        Jff->SetValues(this->row_fsize,this->row_fglobaldoflist,this->row_fsize,this->row_fglobaldoflist,localvalues,ADD_VAL);
    336327
     328                        /*Free ressources:*/
     329                        xfree((void**)&localvalues);
    337330                }
    338                 else{
    339                         _error_(" non dofsymmetrical matrix AddToGlobal routine not support yet!");
    340                 }
    341         #else
    342                 _error_("not supported yet!");
    343         #endif
    344331
     332        }
     333        else{
     334                _error_(" non dofsymmetrical matrix AddToGlobal routine not support yet!");
     335        }
     336
    345337}
    346338/*}}}*/
    347339/*FUNCTION ElementMatrix::CheckConsistency{{{1*/

  • proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Numerics/Vector.cpp

     
    2424/*FUNCTION Vector::Vector(){{{1*/
    2525Vector::Vector(){
    2626
    27         this->M=0;
    2827        #ifdef _HAVE_PETSC_
    2928        this->vector=NULL;
    3029        #else
    … …  
    4140        #ifdef _HAVE_PETSC_
    4241        this->vector=NewVec(pM,fromlocalsize);
    4342        #else
    44         this->M=pM;
    45         this->vector=(double*)xcalloc(pM,sizeof(double));
     43        this->vector=new SeqVec(pM,fromlocalsize);
    4644        #endif
    4745        #ifdef _HAVE_ADOLC_
    4846        this->avector=(adouble*)xmalloc(pM*sizeof(adouble));
    … …  
    5048}
    5149/*}}}*/
    5250/*FUNCTION Vector::Vector(double* serial_vec,int M){{{1*/
    53 Vector::Vector(double* serial_vec,int pM){
     51Vector::Vector(double* serial_vec,int M){
    5452
    5553        int i,j;
    5654
    5755        #ifdef _HAVE_PETSC_
    5856                int* idxm=NULL;
    5957
    60                 this->vector=NewVec(pM);
     58                this->vector=NewVec(M);
    6159               
    6260                idxm=(int*)xmalloc(M*sizeof(int));
    6361                for(i=0;i<M;i++)idxm[i]=i;
    … …  
    6967                xfree((void**)&idxm);
    7068
    7169        #else
    72                 this->M=pM;
    73                 this->vector=(double*)xcalloc(pM,sizeof(double));
     70                this->vector=new SeqVec(serial_vec,M);
    7471        #endif
    7572        #ifdef _HAVE_ADOLC_
    76                 this->avector=(adouble*)xmalloc(pM*sizeof(adouble));
     73                this->avector=(adouble*)xmalloc(M*sizeof(adouble));
    7774        #endif
    7875}
    7976/*}}}*/
    … …  
    8178/*FUNCTION Vector::Vector(Vec petsc_vec){{{1*/
    8279Vector::Vector(Vec petsc_vec){
    8380
    84         /*Get Vector size*/
    85         VecGetSize(petsc_vec,&this->M);
    86 
    8781        /*copy vector*/
    8882        VecDuplicate(petsc_vec,&this->vector);
    8983        VecCopy(petsc_vec,this->vector);
    … …  
    9791        #ifdef _HAVE_PETSC_
    9892        VecFree(&this->vector);
    9993        #else
    100         xfree((void**)&this->vector);
     94        delete this->vector;
    10195        #endif
    10296        #ifdef _HAVE_ADOLC_
    10397        xfree((void**)&this->avector);
    … …  
    114108        #ifdef _HAVE_PETSC_
    115109        VecView(this->vector,PETSC_VIEWER_STDOUT_WORLD);
    116110        #else
    117         printf("Vector size: %i\n",M);
    118         for(i=0;i<M;i++){
    119                 printf("%g\n ",*(vector+i));
    120         }
     111        this->vector->Echo();
    121112        #endif
    122113
    123114        #ifdef _HAVE_ADOLC_
    124         /*Not sure about that one. Should we use the overloaded operator >>?*/
    125         printf("ADOLC Vector equivalent:" );
    126         for(i=0;i<M;i++){
    127                 printf("%g\n ",*(avector+i));
    128         }
     115        /*do nothing for now: */
    129116        #endif
    130117}
    131118/*}}}*/
    … …  
    138125        #ifdef _HAVE_PETSC_
    139126        PetscVectorToMatlabVector(&dataref,this->vector);
    140127        #else
    141         _error_("not implemented yet!");
     128        dataref=this->vector->ToMatlabVector();
    142129        #endif
    143130        return dataref;
    144131
    145132}
    146133/*}}}*/
     134/*FUNCTION MatlabVectorToVector{{{1*/
     135Vector* MatlabVectorToVector(const mxArray* mxvector){
     136
     137        int dummy;
     138        Vector* vector=NULL;
     139
     140        /*allocate vector object: */
     141        vector=new Vector();
     142
     143        #ifdef _HAVE_PETSC_
     144        MatlabVectorToPetscVector(&vector->vector,&dummy,mxvector);
     145        #else
     146        vector->vector=MatlabVectorToSeqVec(mxvector);
     147        #endif
     148       
     149        return vector;
     150}
     151/*}}}*/
    147152#endif
    148153/*FUNCTION Vector::Assemble{{{1*/
    149154void Vector::Assemble(void){
    … …  
    153158                VecAssemblyBegin(this->vector);
    154159                VecAssemblyEnd(this->vector);
    155160        #else
    156                 /*do nothing*/
     161                this->vector->Assemble();
    157162        #endif
    158163
    159164}
    … …  
    166171                _assert_(this->vector);
    167172                VecSetValues(this->vector,ssize,list,values,ISSMToPetscInsertMode(mode));
    168173        #else
    169                 _error_("not implemented yet!");
     174                this->vector->SetValues(ssize,list,values,mode);
    170175        #endif
    171176
    172177}
    … …  
    178183                _assert_(this->vector);
    179184                VecSetValues(this->vector,1,&dof,&value,ISSMToPetscInsertMode(mode));
    180185        #else
    181                 _error_("not implemented yet!");
     186                this->vector->SetValue(dof,value,mode);
    182187        #endif
    183188
    184189}
    … …  
    190195                _assert_(this->vector);
    191196                VecGetValues(this->vector,1,&dof,pvalue);
    192197        #else
    193                 _error_("not implemented yet!");
     198        this->vector->GetValue(pvalue,dof);
    194199        #endif
    195200}
    196201/*}}}*/
    … …  
    201206                _assert_(this->vector);
    202207                VecGetSize(this->vector,pM);
    203208        #else
    204                 _error_("not implemented yet!");
     209                this->vector->GetSize(pM);
    205210        #endif
    206211
    207212}
    … …  
    213218                _assert_(this->vector);
    214219                VecGetLocalSize(this->vector,pM);
    215220        #else
    216                 _error_("not implemented yet!");
     221                this->vector->GetLocalSize(pM);
    217222        #endif
    218223
    219224}
    … …  
    228233                _assert_(this->vector);
    229234                VecDuplicate(this->vector,&vec_output);
    230235                output->vector=vec_output;
    231                 VecGetSize(output->vector,&output->M);
    232236        #else
    233                 _error_("not implemented yet!");
     237                output->vector=this->vector->Duplicate();
    234238        #endif
    235239        return output;
    236240
    … …  
    243247                _assert_(this->vector);
    244248                VecSet(this->vector,value);
    245249        #else
    246                 _error_("not implemented yet!");
     250                this->vector->Set(value);
    247251        #endif
    248252
    249253}
    … …  
    255259                _assert_(this->vector);
    256260                VecAXPY(this->vector,a,X->vector);
    257261        #else
    258                 _error_("not implemented yet!");
     262                this->vector->AXPY(X->vector,a);
    259263        #endif
    260264}
    261265/*}}}*/
    … …  
    266270                _assert_(this->vector);
    267271                VecAYPX(this->vector,a,X->vector);
    268272        #else
    269                 _error_("not implemented yet!");
     273                this->vector->AYPX(X->vector,a);
    270274        #endif
    271275
    272276}
    … …  
    279283        #ifdef _HAVE_PETSC_
    280284                VecToMPISerial(&vec_serial, this->vector);
    281285        #else
    282                 _error_("not implemented yet!");
     286                vec_serial=this->vector->ToMPISerial();
    283287        #endif
    284288
    285289        return vec_serial;
    … …  
    292296        #ifdef _HAVE_PETSC_
    293297                VecCopy(this->vector,to->vector);
    294298        #else
    295                 _error_("not implemented yet!");
     299                this->vector->Copy(to->vector);
    296300        #endif
    297301
    298302}
    … …  
    305309                _assert_(this->vector);
    306310                VecNorm(this->vector,ISSMToPetscNormMode(norm_type),&norm);
    307311        #else
    308                 _error_("not implemented yet!");
     312                norm=this->vector->Norm(norm_type);
    309313        #endif
    310314        return norm;
    311315}
    … …  
    317321                _assert_(this->vector);
    318322                VecScale(this->vector,scale_factor);
    319323        #else
    320                 _error_("not implemented yet!");
     324                this->vector->Scale(scale_factor);
    321325        #endif
    322326}
    323327/*}}}*/
    … …  
    329333                _assert_(this->vector);
    330334                VecDot(this->vector,vector->vector,&dot);
    331335        #else
    332                 _error_("not implemented yet!");
     336                dot=this->vector->Dot(vector->vector);
    333337        #endif
    334338        return dot;
    335339}
    … …  
    341345                _assert_(this->vector);
    342346                VecPointwiseDivide(this->vector,x->vector,y->vector);
    343347        #else
    344                 _error_("not implemented yet!");
     348                this->vector->PointwiseDivide(x->vector,y->vector);
    345349        #endif
    346350}
    347351/*}}}*/

  • proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Numerics/Matrix.cpp

     
    4747        #ifdef _HAVE_PETSC_
    4848        this->matrix=NewMat(pM,pN);
    4949        #else
    50         this->matrix=(double*)xcalloc(pM*pN,sizeof(double));
     50        this->matrix=new SeqMat(pM,pN);
    5151        #endif
    5252        #ifdef _HAVE_ADOLC_
    5353        this->amatrix=(adouble*)xmalloc(pM*pN*sizeof(adouble));
    … …  
    6363        #ifdef _HAVE_PETSC_
    6464        this->matrix=NewMat(pM,pN,sparsity);
    6565        #else
    66         this->matrix=(double*)xcalloc(pM*pN,sizeof(double));
     66        this->matrix=new SeqMat(pM,pN,sparsity);
    6767        #endif
    6868        #ifdef _HAVE_ADOLC_
    6969        this->amatrix=(adouble*)xmalloc(pM*pN*sizeof(adouble));
    … …  
    9696        xfree((void**)&idxm);
    9797        xfree((void**)&idxn);
    9898        #else
    99         this->matrix=(double*)xcalloc(pM*pN,sizeof(double));
     99        this->matrix=new SeqMat(serial_mat,pM,pN,sparsity);
    100100        #endif
    101101        #ifdef _HAVE_ADOLC_
    102102        this->amatrix=(adouble*)xmalloc(pM*pN*sizeof(adouble));
    … …  
    112112        #ifdef _HAVE_PETSC_
    113113        this->matrix=NewMat(pM,pN,connectivity,numberofdofspernode);
    114114        #else
    115         this->matrix=(double*)xcalloc(pM*pN,sizeof(double));
     115        this->matrix=new SeqMat(pM,pN,connectivity,numberofdofspernode);
    116116        #endif
    117117        #ifdef _HAVE_ADOLC_
    118118        this->amatrix=(adouble*)xmalloc(pM*pN*sizeof(adouble));
    … …  
    125125        #ifdef _HAVE_PETSC_
    126126        MatFree(&this->matrix);
    127127        #else
    128         xfree((void**)&this->matrix);
     128        delete this->matrix;
    129129        #endif
    130130        #ifdef _HAVE_ADOLC_
    131131        xfree((void**)&this->amatrix);
    … …  
    142142        #ifdef _HAVE_PETSC_
    143143        MatView(this->matrix,PETSC_VIEWER_STDOUT_WORLD);
    144144        #else
    145         printf("Matrix size: %i-%i\n",M,N);
    146         for(i=0;i<M;i++){
    147                 for(j=0;j<N;j++){
    148                         printf("%g ",*(matrix+N*i+j));
    149                 }
    150                 printf("\n");
    151         }
     145        this->matrix->Echo();
    152146        #endif
    153147
    154148        #ifdef _HAVE_ADOLC_
    … …  
    172166        #ifdef _HAVE_PETSC_
    173167        PetscMatrixToMatlabMatrix(&dataref,this->matrix);
    174168        #else
    175         _error_("not implemented yet!");
     169        dataref=this->matrix->ToMatlabMatrix();
    176170        #endif
    177171        return dataref;
    178172
    179173}
    180174/*}}}*/
     175/*FUNCTION MatlabMatrixToMatrix{{{1*/
     176Matrix* MatlabMatrixToMatrix(const mxArray* mxmatrix){
     177
     178        int dummy;
     179        Matrix* matrix=NULL;
     180
     181        /*allocate matrix object: */
     182        matrix=new Matrix();
     183
     184        #ifdef _HAVE_PETSC_
     185        MatlabMatrixToPetscMatrix(&matrix->matrix,&matrix->M,&matrix->N,mxmatrix);
     186        #else
     187        matrix->matrix=MatlabMatrixToSeqMat(mxmatrix);
     188        #endif
     189       
     190        return matrix;
     191}
     192/*}}}*/
    181193#endif
    182194/*FUNCTION Matrix::Assemble{{{1*/
    183195void Matrix::Assemble(void){
    … …  
    189201                        MatCompress(this->matrix);
    190202                #endif
    191203        #else
    192                 /*do nothing:*/
     204                this->matrix->Assemble();
    193205        #endif
    194206
    195207}
    … …  
    202214                _assert_(this->matrix);
    203215                MatNorm(this->matrix,ISSMToPetscNormMode(norm_type),&norm);
    204216        #else
    205                 _error_("not implemented yet!");
     217                norm=this->matrix->Norm(norm_type);
    206218        #endif
    207219        return norm;
    208220}
    … …  
    214226                _assert_(this->matrix);
    215227                MatGetSize(this->matrix,pM,pN);
    216228        #else
    217                 _error_("not implemented yet!");
     229                this->matrix->GetSize(pM,pN);
    218230        #endif
    219231}
    220232/*}}}*/
    … …  
    225237                _assert_(this->matrix);
    226238                MatGetLocalSize(this->matrix,pM,pN);
    227239        #else
    228                 _error_("not implemented yet!");
     240                this->matrix->GetLocalSize(pM,pN);
    229241        #endif
    230242}
    231243/*}}}*/
    … …  
    237249                _assert_(X->vector);
    238250                MatMultPatch(this->matrix,X->vector,AX->vector);
    239251        #else
    240                 _error_("not implemented yet!");
     252                this->matrix->MatMult(X->vector,AX->vector);
    241253        #endif
    242254}
    243255/*}}}*/
    … …  
    248260
    249261        output=new Matrix();
    250262
    251         #ifdef _HAVE_PETSC
     263        #ifdef _HAVE_PETSC_
    252264                _assert_(this->matrix);
    253265                MatDuplicate(this->matrix,MAT_COPY_VALUES,&output->matrix);
    254266        #else
    255                 _error_("not implemented yet!");
     267                output->matrix=this->matrix->Duplicate();
    256268        #endif
    257269}
    258270/*}}}*/
    … …  
    261273
    262274        double* output=NULL;
    263275
    264         #ifdef _HAVE_PETSC
     276        #ifdef _HAVE_PETSC_
    265277                MatToSerial(&output,this->matrix);
    266278        #else
    267                 _error_("not implemented yet!");
     279                output=this->matrix->ToSerial();
    268280        #endif
    269281        return output;
    270282}
    271283/*}}}*/
     284/*FUNCTION Matrix::SetValues{{{1*/
     285void Matrix::SetValues(int m,int* idxm,int n,int* idxn,double* values,InsMode mode){
     286
     287        #ifdef _HAVE_PETSC_
     288                MatSetValues(this->matrix,m,idxm,n,idxn,values,ISSMToPetscInsertMode(mode));
     289        #else
     290                this->matrix->SetValues(m,idxm,n,idxn,values,mode);
     291        #endif
     292}
     293/*}}}*/

  • proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Numerics/Vector.h

     
    3030
    3131        public:
    3232       
    33                 int M;
    34 
    3533                #ifdef _HAVE_PETSC_
    3634                Vec vector;
    3735                #else
    38                 double* vector;
     36                SeqVec* vector;
    3937                #endif
    4038                #ifdef _HAVE_ADOLC_
    4139                adouble* avector;
    … …  
    7371                double Dot(Vector* vector);
    7472                /*}}}*/
    7573};
     74
     75/*API: */
     76#ifdef _SERIAL_
     77Vector* MatlabVectorToVector(const mxArray* mxvector);
     78#endif
     79
    7680#endif //#ifndef _VECTOR_H_

  • proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Numerics/Matrix.h

     
    3636                #ifdef _HAVE_PETSC_
    3737                Mat matrix;
    3838                #else
    39                 double* matrix;
     39                SeqMat* matrix;
    4040                #endif
    4141                #ifdef _HAVE_ADOLC_
    4242                adouble* amatrix;
    … …  
    6262                void MatMult(Vector* X,Vector* AX);
    6363                Matrix* Duplicate(void);
    6464                double* ToSerial(void);
     65                void SetValues(int m,int* idxm,int n,int* idxn,double* values,InsMode mode);
    6566                /*}}}*/
    6667
    6768};
     69/*API: */
     70#ifdef _SERIAL_
     71Matrix* MatlabMatrixToMatrix(const mxArray* mxmatrix);
     72#endif
     73
    6874#endif //#ifndef _MATRIX_H_

  • proj/ice/larour/issm-uci-clean/trunk-jpl/src/c/objects/Numerics/ElementVector.cpp

     
    165165        int i;
    166166        double* localvalues=NULL;
    167167
    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);
     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]];
    180173                }
    181         #else
    182                 _error_("not supported yet!");
    183         #endif
     174                /*add local values into global  vector, using the fglobaldoflist: */
     175                pf->SetValues(this->fsize,this->fglobaldoflist,localvalues,ADD_VAL);
     176
     177                /*Free ressources:*/
     178                xfree((void**)&localvalues);
     179        }
    184180       
    185181}
    186182/*}}}*/
    … …  
    190186        int i;
    191187        double* localvalues=NULL;
    192188
    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);
     189        if(this->fsize){
     190                /*first, retrieve values that are in the f-set from the g-set values vector: */
     191                localvalues=(double*)xmalloc(this->fsize*sizeof(double));
     192                for(i=0;i<this->fsize;i++){
     193                        localvalues[i]=this->values[this->flocaldoflist[i]];
    205194                }
    206         #else
    207                 _error_("not supported yet!");
    208         #endif
     195                /*add local values into global  vector, using the fglobaldoflist: */
     196                pf->SetValues(this->fsize,this->fglobaldoflist,localvalues,INS_VAL);
    209197
     198                /*Free ressources:*/
     199                xfree((void**)&localvalues);
     200        }
     201
    210202}
    211203/*}}}*/
    212204/*FUNCTION ElementVector::Echo{{{1*/
Note: See TracBrowser for help on using the repository browser.