Changeset 25752


Ignore:
Timestamp:
11/13/20 16:20:12 (4 years ago)
Author:
Eric.Larour
Message:

CHG: fixed issues with barystatic partitioning that was missing from the prototypes.

Location:
issm/trunk-jpl/src/c
Files:
9 edited

Legend:

Unmodified
Added
Removed

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

    r25627 r25752  
    377377                virtual IssmDouble    OceanAverage(IssmDouble* Sg, SealevelMasks* masks)=0;
    378378                virtual void          SealevelriseMomentOfInertia(IssmDouble* dI_list,IssmDouble* Sg_old, SealevelMasks* masks)=0;
    379                 virtual IssmDouble    SealevelriseEustaticIce(IssmDouble* Sgi, SealevelMasks* masks,IssmDouble oceanarea)=0;
    380                 virtual IssmDouble    SealevelriseEustaticHydro(IssmDouble* Sgi, SealevelMasks* masks,IssmDouble oceanarea)=0;
     379                virtual IssmDouble    SealevelriseEustaticIce(IssmDouble* Sgi, SealevelMasks* masks,Vector<IssmDouble>* barystatic_contribution,IssmDouble* partition,IssmDouble oceanarea)=0;
     380                virtual IssmDouble    SealevelriseEustaticHydro(IssmDouble* Sgi, SealevelMasks* masks,Vector<IssmDouble>* barystatic_contribution, IssmDouble* partition, IssmDouble oceanarea)=0;
    381381                virtual void          SealevelriseEustaticBottomPressure(IssmDouble* Sgi, SealevelMasks* masks)=0;
    382382                virtual void          SealevelriseGeometry(IssmDouble* latitude,IssmDouble* longitude,IssmDouble* radius,IssmDouble* xx, IssmDouble* yy, IssmDouble* zz)=0;

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

    r25627 r25752  
    217217                void    SealevelriseMomentOfInertia(IssmDouble* dI_list,IssmDouble* Sg_old, SealevelMasks* masks){_error_("not implemented yet!");};
    218218                void    SealevelriseGeometry(IssmDouble* latitude,IssmDouble* longitude,IssmDouble* radius,IssmDouble* xx, IssmDouble* yy, IssmDouble* zz){_error_("not implemented yet!");};
    219                 IssmDouble    SealevelriseEustaticIce(IssmDouble* Sgi, SealevelMasks* masks, IssmDouble oceanarea){_error_("not implemented yet!");};
    220                 IssmDouble    SealevelriseEustaticHydro(IssmDouble* Sgi, SealevelMasks* masks, IssmDouble oceanarea){_error_("not implemented yet!");};
     219                IssmDouble    SealevelriseEustaticIce(IssmDouble* Sgi, SealevelMasks* masks, Vector<IssmDouble>* barystatic_contribution, IssmDouble* partition, IssmDouble oceanarea){_error_("not implemented yet!");};
     220                IssmDouble    SealevelriseEustaticHydro(IssmDouble* Sgi, SealevelMasks* masks, Vector<IssmDouble>* barystatic_contribution, IssmDouble* partition,IssmDouble oceanarea){_error_("not implemented yet!");};
    221221                void    SealevelriseEustaticBottomPressure(IssmDouble* Sgi,SealevelMasks* masks){_error_("not implemented yet!");};
    222222                void    SealevelriseNonEustatic(IssmDouble* Sgo, IssmDouble* Sg_old, SealevelMasks* masks){_error_("not implemented yet!");};

  • issm/trunk-jpl/src/c/classes/Elements/Seg.h

    r25627 r25752  
    173173                void    SetSealevelMasks(SealevelMasks* masks){_error_("not implemented yet!");};
    174174                void    SealevelriseGeometry(IssmDouble* latitude,IssmDouble* longitude,IssmDouble* radius,IssmDouble* xx, IssmDouble* yy, IssmDouble* zz){_error_("not implemented yet!");};
    175                 IssmDouble    SealevelriseEustaticIce(IssmDouble* Sgi, SealevelMasks* masks, IssmDouble oceanarea){_error_("not implemented yet!");};
    176                 IssmDouble    SealevelriseEustaticHydro(IssmDouble* Sgi, SealevelMasks* masks, IssmDouble oceanarea){_error_("not implemented yet!");};
     175                IssmDouble    SealevelriseEustaticIce(IssmDouble* Sgi, SealevelMasks* masks, Vector<IssmDouble>* barystatic_contribution, IssmDouble* partition, IssmDouble oceanarea){_error_("not implemented yet!");};
     176                IssmDouble    SealevelriseEustaticHydro(IssmDouble* Sgi, SealevelMasks* masks, Vector<IssmDouble>* barystatic_contribution, IssmDouble* partition,IssmDouble oceanarea){_error_("not implemented yet!");};
    177177                void    SealevelriseEustaticBottomPressure(IssmDouble* Sgi,SealevelMasks* masks){_error_("not implemented yet!");};
    178178                void    SealevelriseNonEustatic(IssmDouble* Sgo,IssmDouble* Sg_old, SealevelMasks* masks){_error_("not implemented yet!");};

  • issm/trunk-jpl/src/c/classes/Elements/Tetra.h

    r25627 r25752  
    179179                void    SealevelriseMomentOfInertia(IssmDouble* dI_list,IssmDouble* Sg_old, SealevelMasks* masks){_error_("not implemented yet!");};
    180180                void    SealevelriseGeometry(IssmDouble* latitude,IssmDouble* longitude,IssmDouble* radius,IssmDouble* xx, IssmDouble* yy, IssmDouble* zz){_error_("not implemented yet!");};
    181                 IssmDouble    SealevelriseEustaticIce(IssmDouble* Sgi, SealevelMasks* masks, IssmDouble oceanarea){_error_("not implemented yet!");};
    182                 IssmDouble    SealevelriseEustaticHydro(IssmDouble* Sgi, SealevelMasks* masks, IssmDouble oceanarea){_error_("not implemented yet!");};
     181                IssmDouble    SealevelriseEustaticIce(IssmDouble* Sgi, SealevelMasks* masks, Vector<IssmDouble>* barystatic_contribution, IssmDouble* partition, IssmDouble oceanarea){_error_("not implemented yet!");};
     182                IssmDouble    SealevelriseEustaticHydro(IssmDouble* Sgi, SealevelMasks* masks, Vector<IssmDouble>* barystatic_contribution, IssmDouble* partition,IssmDouble oceanarea){_error_("not implemented yet!");};
    183183                void    SealevelriseEustaticBottomPressure(IssmDouble* Sgi,SealevelMasks* masks){_error_("not implemented yet!");};
    184184                void    SealevelriseNonEustatic(IssmDouble* Sgo, IssmDouble* Sg_old, SealevelMasks* masks){_error_("not implemented yet!");};

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

    r25751 r25752  
    55205520                /*Rigid earth gravitational perturbation: */
    55215521                G[i]+=G_rigid_precomputed[index];
    5522                 if(elastic) G[i]+=G_elastic_precomputed[index];
     5522                if(computeelastic) G[i]+=G_elastic_precomputed[index];
    55235523                G[i]=G[i]*constant;
    55245524
    … …  
    55705570        #else
    55715571        xDelete(G);
    5572         if(elastic){
     5572        if(computeelastic){
    55735573                xDelete(GU);
    55745574                if(horiz){
    … …  
    55825582}
    55835583/*}}}*/
    5584 IssmDouble Tria::SealevelriseEustaticIce(IssmDouble* Sgi, SealevelMasks* masks, IssmDouble oceanarea){ /*{{{*/
     5584IssmDouble Tria::SealevelriseEustaticIce(IssmDouble* Sgi, SealevelMasks* masks, Vector<IssmDouble>* barystatic_contribution, IssmDouble* partition, IssmDouble oceanarea){ /*{{{*/
    55855585
    55865586        /*diverse:*/
    … …  
    56545654                if(glfraction==0)phi=1;
    56555655                #ifdef _ISSM_DEBUG_
    5656                 this->AddInput2(SealevelEustaticMaskEnum,&phi,P0Enum);
     5656                this->AddInput(SealevelEustaticMaskEnum,&phi,P0Enum);
    56575657                #endif
    56585658        }
    … …  
    57125712}
    57135713/*}}}*/
    5714 IssmDouble Tria::SealevelriseEustaticHydro(IssmDouble* Sgi, SealevelMasks* masks, IssmDouble oceanarea){ /*{{{*/
     5714IssmDouble Tria::SealevelriseEustaticHydro(IssmDouble* Sgi, SealevelMasks* masks, Vector<IssmDouble>* barystatic_contribution, IssmDouble* partition, IssmDouble oceanarea){ /*{{{*/
    57155715
    57165716        /*diverse:*/

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

    r25627 r25752  
    165165                void    SetSealevelMasks(SealevelMasks* masks);
    166166                void    SealevelriseGeometry(IssmDouble* latitude,IssmDouble* longitude,IssmDouble* radius, IssmDouble* xx, IssmDouble* yy, IssmDouble* zz);
    167                 IssmDouble    SealevelriseEustaticIce(IssmDouble* Sgi, SealevelMasks* masks, IssmDouble oceanarea);
    168                 IssmDouble    SealevelriseEustaticHydro(IssmDouble* Sgi, SealevelMasks* masks, IssmDouble oceanarea);
     167                IssmDouble    SealevelriseEustaticIce(IssmDouble* Sgi, SealevelMasks* masks, Vector<IssmDouble>* barystatic_contribution, IssmDouble* partition, IssmDouble oceanarea);
     168                IssmDouble    SealevelriseEustaticHydro(IssmDouble* Sgi, SealevelMasks* masks, Vector<IssmDouble>* barystatic_contribution, IssmDouble* partition,IssmDouble oceanarea);
    169169                void    SealevelriseEustaticBottomPressure(IssmDouble* Sgi,SealevelMasks* masks);
    170170                void    SealevelriseNonEustatic(IssmDouble* Sgo,IssmDouble* Sg_old,SealevelMasks* masks);

  • issm/trunk-jpl/src/c/classes/FemModel.cpp

    r25727 r25752  
    47484748#endif
    47494749#ifdef _HAVE_SEALEVELRISE_
    4750 void FemModel::SealevelriseEustatic(Vector<IssmDouble>* pRSLgi, IssmDouble* poceanarea, IssmDouble* pbslr,IssmDouble* pbslrice, IssmDouble* pbslrhydro, SealevelMasks* masks) { /*{{{*/
     4750void FemModel::SealevelriseEustatic(Vector<IssmDouble>* pRSLgi, IssmDouble* poceanarea, IssmDouble* pbslr,IssmDouble* pbslrice, IssmDouble* pbslrhydro, IssmDouble** pbslrice_partition,IssmDouble** pbslrhydro_partition,SealevelMasks* masks) { /*{{{*/
    47514751
    47524752        /*serialized vectors:*/
    … …  
    47604760        int bp_compute_fingerprints= 0;
    47614761
     4762        Vector<IssmDouble>* bslrice_partition=NULL;
     4763        IssmDouble* bslrice_partition_serial=NULL;
     4764        IssmDouble* partitionice=NULL;
     4765        int npartice,nel;
     4766
     4767        Vector<IssmDouble>* bslrhydro_partition=NULL;
     4768        IssmDouble* bslrhydro_partition_serial=NULL;
     4769        IssmDouble* partitionhydro=NULL;
     4770        int nparthydro;
     4771               
     4772
    47624773   /*Initialize temporary vector that will be used to sum eustatic components
    47634774    * on all local elements, prior to assembly:*/
    … …  
    47794790        _assert_(oceanarea>0.);
    47804791
     4792        /*Initialize partition vectors to retrieve barystatic contributions: */
     4793        this->parameters->FindParam(&npartice,SolidearthNpartIceEnum);
     4794        if(npartice){
     4795                this->parameters->FindParam(&partitionice,&nel,NULL,SolidearthPartitionIceEnum);
     4796                bslrice_partition= new Vector<IssmDouble>(npartice);
     4797        }
     4798
     4799        this->parameters->FindParam(&nparthydro,SolidearthNpartHydroEnum);
     4800        if(nparthydro){
     4801                this->parameters->FindParam(&partitionhydro,&nel,NULL,SolidearthPartitionHydroEnum);
     4802                bslrhydro_partition= new Vector<IssmDouble>(nparthydro);
     4803        }
     4804
     4805
    47814806        /*Call the sea level rise core for ice : */
    47824807        bslrice_cpu=0;
    47834808        for(Object* & object : this->elements->objects){
    47844809                Element* element = xDynamicCast<Element*>(object);
    4785                 bslrice_cpu+=element->SealevelriseEustaticIce(RSLgi,masks, oceanarea);
     4810                bslrice_cpu+=element->SealevelriseEustaticIce(RSLgi,masks, bslrice_partition,partitionice,oceanarea);
    47864811        }
    47874812
    … …  
    47904815        for(int i=0;i<elements->Size();i++){
    47914816                Element* element=xDynamicCast<Element*>(elements->GetObjectByOffset(i));
    4792                 bslrhydro_cpu+=element->SealevelriseEustaticHydro(RSLgi,masks, oceanarea);
     4817                bslrhydro_cpu+=element->SealevelriseEustaticHydro(RSLgi,masks, bslrhydro_partition,partitionhydro,oceanarea);
    47934818        }
    47944819
    … …  
    48154840        _assert_(!xIsNan<IssmDouble>(bslrhydro));
    48164841
     4842        /*Take care of partition vectors:*/
     4843        if(bslrice_partition){
     4844                bslrice_partition->Assemble();
     4845                bslrice_partition_serial=bslrice_partition->ToMPISerial();
     4846        }
     4847        if(bslrhydro_partition){
     4848                bslrhydro_partition->Assemble();
     4849                bslrhydro_partition_serial=bslrhydro_partition->ToMPISerial();
     4850        }
     4851
     4852
    48174853        /*Free ressources:*/
    48184854        xDelete<int>(indices);
    48194855        xDelete<IssmDouble>(RSLgi);
     4856        if(bslrhydro_partition)delete bslrhydro_partition;
     4857        if(bslrice_partition)delete bslrice_partition;
     4858        if(partitionhydro)xDelete<IssmDouble>(partitionhydro);
     4859        if(partitionice)xDelete<IssmDouble>(partitionice);
    48204860
    48214861        /*Assign output pointers:*/
    … …  
    48244864        *pbslrhydro  = bslrhydro;
    48254865        *pbslr=bslrhydro+bslrice;
     4866        *pbslrice_partition=bslrice_partition_serial;
     4867        *pbslrhydro_partition=bslrhydro_partition_serial;
    48264868
    48274869}
    … …  
    48374879
    48384880        bool computerigid = true;
     4881        bool computeelastic= true;
    48394882
    48404883        /*recover computational flags: */

  • issm/trunk-jpl/src/c/classes/FemModel.h

    r25627 r25752  
    166166                #endif
    167167                #ifdef _HAVE_SEALEVELRISE_
    168                 void SealevelriseEustatic(Vector<IssmDouble>* pRSLgi, IssmDouble* poceanarea, IssmDouble* pbslr,IssmDouble* pbslrice, IssmDouble* pbslrhydro, SealevelMasks* masks);
     168                void SealevelriseEustatic(Vector<IssmDouble>* pRSLgi, IssmDouble* poceanarea, IssmDouble* pbslr,IssmDouble* pbslrice, IssmDouble* pbslrhydro, IssmDouble** pbslrice_partition,IssmDouble** pbslrhydro_partition, SealevelMasks* masks);
    169169                void SealevelriseNonEustatic(Vector<IssmDouble>* pSgo, Vector<IssmDouble>* pSg_old,  SealevelMasks* masks,bool verboseconvolution);
    170170                void SealevelriseRotationalFeedback(Vector<IssmDouble>* pRSLgo_rot, Vector<IssmDouble>* pRSLg_old, IssmDouble* pIxz, IssmDouble* pIyz, IssmDouble* pIzz, SealevelMasks* masks);

  • issm/trunk-jpl/src/c/cores/sealevelchange_core.cpp

    r25680 r25752  
    4040        if(!isslr)return;
    4141
     42        /*Verbose: */
     43        if(VerboseSolution()) _printf0_("   computing sea level rise\n");
     44
    4245        /*Run gia core: */
    4346        if(isgia){
    … …  
    6568                int     numoutputs;
    6669                char **requested_outputs = NULL;
     70                if(VerboseSolution()) _printf0_("   saving results\n");
    6771                femmodel->parameters->FindParam(&requested_outputs,&numoutputs,SealevelriseRequestedOutputsEnum);
    6872                femmodel->RequestedOutputsx(&femmodel->results,requested_outputs,numoutputs);
    … …  
    430434        IssmDouble bslr;
    431435        IssmDouble bslrice;
     436        IssmDouble* bslrice_partition=NULL;
    432437        IssmDouble bslrhydro;
     438        IssmDouble* bslrhydro_partition=NULL;
    433439        IssmDouble cumbslr;
    434440        IssmDouble cumbslrice;
    435441        IssmDouble cumbslrhydro;
     442        IssmDouble* cumbslrice_partition=NULL;
     443        int npartice;
     444        IssmDouble* cumbslrhydro_partition=NULL;
     445        int nparthydro;
    436446       
    437447        if(VerboseSolution()) _printf0_("         computing bslr components on ice\n");
    … …  
    446456        femmodel->parameters->FindParam(&cumbslrice,CumBslrIceEnum);
    447457        femmodel->parameters->FindParam(&cumbslrhydro,CumBslrHydroEnum);
     458        femmodel->parameters->FindParam(&npartice,SolidearthNpartIceEnum);
     459        femmodel->parameters->FindParam(&nparthydro,SolidearthNpartHydroEnum);
     460        if(npartice) femmodel->parameters->FindParam(&cumbslrice_partition,&npartice,NULL,CumBslrIcePartitionEnum);
     461        if(nparthydro) femmodel->parameters->FindParam(&cumbslrhydro_partition,&nparthydro,NULL,CumBslrHydroPartitionEnum);
    448462
    449463        /*Initialize:*/
    … …  
    451465
    452466        /*call the bslr main module: */
    453         femmodel->SealevelriseEustatic(RSLgi,&oceanarea,&bslr, &bslrice, &bslrhydro, masks); //this computes
     467        femmodel->SealevelriseEustatic(RSLgi,&oceanarea,&bslr, &bslrice, &bslrhydro, &bslrice_partition, &bslrhydro_partition,masks); //this computes
    454468
    455469        /*we need to average RSLgi over the ocean: RHS term  4 in Eq.4 of Farrel and clarke. Only the elements can do that: */
    … …  
    465479        femmodel->results->AddResult(new GenericExternalResult<IssmDouble>(femmodel->results->Size()+1,BslrHydroEnum,-bslrhydro,step,time));
    466480
     481        //cumulative barystatic contribution:
    467482        cumbslr=cumbslr-bslr;
    468483        cumbslrice=cumbslrice-bslrice;
    469484        cumbslrhydro=cumbslrhydro-bslrhydro;
     485
    470486        femmodel->parameters->SetParam(cumbslr,CumBslrEnum);
    471487        femmodel->parameters->SetParam(cumbslrice,CumBslrIceEnum);
    … …  
    474490        femmodel->results->AddResult(new GenericExternalResult<IssmDouble>(femmodel->results->Size()+1,CumBslrEnum,cumbslr,step,time));
    475491        femmodel->results->AddResult(new GenericExternalResult<IssmDouble>(femmodel->results->Size()+1,CumBslrIceEnum,cumbslrice,step,time));
    476         femmodel->results->AddResult(new GenericExternalResult<IssmDouble>(femmodel->results->Size()+1,CumBslrHydroEnum,cumbslrhydro,step,time));
     492        femmodel->results->AddResult(new GenericExternalResult<IssmDouble>(femmodel->results->Size()+1,CumBslrHydroEnum,cumbslrhydro,step,time))
     493
     494        //cumulative barystatic contributions by partition:
     495        if(npartice){
     496                for(int i=0;i<npartice;i++) cumbslrice_partition[i] -= bslrice_partition[i];
     497                femmodel->parameters->SetParam(cumbslrice_partition,npartice,1,CumBslrIcePartitionEnum);
     498                femmodel->results->AddResult(new GenericExternalResult<IssmPDouble*>(femmodel->results->Size()+1,CumBslrIcePartitionEnum,cumbslrice_partition,npartice,1,step,time));
     499        }
     500
     501        if(nparthydro){
     502                for(int i=0;i<nparthydro;i++) cumbslrhydro_partition[i] -= bslrhydro_partition[i];
     503                femmodel->parameters->SetParam(cumbslrhydro_partition,nparthydro,1,CumBslrHydroPartitionEnum);
     504                femmodel->results->AddResult(new GenericExternalResult<IssmPDouble*>(femmodel->results->Size()+1,CumBslrHydroPartitionEnum,cumbslrhydro_partition,nparthydro,1,step,time));
     505        }
    477506        /*}}}*/
    478507       
Note: See TracChangeset for help on using the changeset viewer.