Changeset 8549


Ignore:
Timestamp:
06/08/11 11:57:29 (14 years ago)
Author:
bkhakbaz
Message:

added the hydr_CR and hydro_n as parameters of hydrology and remove hydro_gamma from hydrology parameters

Location:
issm/trunk/src
Files:
9 edited

Legend:

Unmodified
Added
Removed

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

    r8540 r8549  
    991991        double     Ke_gg_gaussian2[numdof][numdof]  ={0.0};
    992992        double     Ke_gg_gaussian3[numdof][numdof]  ={0.0};
    993         double     hydro_gamma,hydro_p,hydro_q;
     993        double     hydro_gamma,hydro_p,hydro_q,hydro_CR,hydro_n,rho_water,g;
    994994        double     DL_scalar1;
    995995        double     DL_scalar2;
    … …  
    10051005
    10061006        /*retrieve material parameters: */
    1007         hydro_gamma=matpar->GetGamma();
     1007//      hydro_gamma=matpar->GetGamma();
     1008        g          =matpar->GetG();
     1009        rho_water  =matpar->GetRhoWater();
    10081010        hydro_p    =matpar->GetHydroP();
    10091011        hydro_q    =matpar->GetHydroQ();
    1010 
     1012   hydro_CR    =matpar->GetHydroCR();
     1013        hydro_n    =matpar->GetHydroN();
     1014       
     1015        hydro_gamma=1/(hydro_n*pow(hydro_CR,hydro_p)*pow((rho_water*g),hydro_q));
     1016       
    10111017        /* Start  looping on the number of gaussian points: */
    10121018        gauss=new GaussTria(2);
    … …  
    10181024                GetNodalFunctions(&L[0],gauss);
    10191025                GetNodalFunctionsDerivatives(&dL[0][0],&xyz_list[0][0],gauss);
    1020 
     1026 //printf("NUMVERTICES=%d, Jdettria=%f\n",NUMVERTICES,Jdettria); //BK
    10211027                for(i=0;i<NUMVERTICES;i++)dLx[i]=dL[0][i];
    1022                 for(i=0;i<NUMVERTICES;i++)dLy[i]=dL[1][i];
     1028           for(i=0;i<NUMVERTICES;i++)dLy[i]=dL[1][i];
     1029  /*    for(i=0;i<NUMVERTICES;i++)
     1030                  { dLx[i]=dL[0][i];
     1031                printf("dLx[%d]=%f, L[%d]=%f\n",i,dLx[i]),i,L[i];}
     1032                for(i=0;i<NUMVERTICES;i++)
     1033                         { dLy[i]=dL[1][i];
     1034                                               printf("dLy[%d]=%f\n",i,dLy[i]);} */
    10231035
    10241036                /*Get K parameter: */
    10251037                GetHydrologyK(&K[0],&xyz_list[0][0],gauss);
    1026        
     1038//      printf("K[0]=%g,K[1]=%g\n",K[0],K[1]);
     1039
    10271040                if(dt){
    10281041                        DL_scalar1=gauss->weight*Jdettria;
    … …  
    10491062                for(i=0;i<numdof;i++) for(j=0;j<numdof;j++) Ke->values[i*numdof+j]+=Ke_gg_gaussian2[i][j];
    10501063                for(i=0;i<numdof;i++) for(j=0;j<numdof;j++) Ke->values[i*numdof+j]+=Ke_gg_gaussian3[i][j];
     1064
    10511065        }
    10521066
    … …  
    20912105                basal_melting_input->GetParameterValue(&basal_melting_g,gauss);
    20922106                old_watercolumn_input->GetParameterValue(&old_watercolumn_g,gauss);
    2093        
     2107//      printf("ig=%d, old_watercolumn_g=%f,Jdettria=%f\n",ig,old_watercolumn_g,Jdettria);//bk
     2108//      if(old_watercolumn_g<0) printf("old_watercolumn_g=%g,Jdettria=%f\n",old_watercolumn_g,Jdettria);//bk
     2109
    20942110                if(dt)for(i=0;i<numdof;i++) pe->values[i]+=Jdettria*gauss->weight*(old_watercolumn_g+dt*basal_melting_g)*L[i];
    20952111                else  for(i=0;i<numdof;i++) pe->values[i]+=Jdettria*gauss->weight*basal_melting_g*L[i];
    … …  
    26992715                watercolumn_input->GetParameterValue(&watercolumn,gauss);
    27002716                values[i]=watercolumn;
    2701                 if((values[i])<0)values[i]=0;
     2717        //      if((values[i])<0)values[i]=0;
    27022718        }
    27032719
    … …  
    27852801        K[0]=pow(w,2)*(rho_ice*g*dsdx+(rho_water/rho_ice-1)*rho_ice*g*dbdx) - rho_ice * g * kn* w * (dsdx - dbdx ) * surface_slope;
    27862802        K[1]=pow(w,2)*(rho_ice*g*dsdy+(rho_water/rho_ice-1)*rho_ice*g*dbdy) - rho_ice * g * kn *w * (dsdy - dbdy ) * surface_slope;
    2787 
    2788 }
     2803 
     2804        //bk
     2805          // K[0]=fabs(pow(w,2)*(rho_ice*g*dsdx+(rho_water/rho_ice-1)*rho_ice*g*dbdx) - rho_ice * g * kn* w * (dsdx - dbdx ) * surface_slope);
     2806          // K[1]=fabs(pow(w,2)*(rho_ice*g*dsdy+(rho_water/rho_ice-1)*rho_ice*g*dbdy) - rho_ice * g * kn *w * (dsdy - dbdy ) * surface_slope);
     2807//       printf("K[0]=%g,K[1]=%g,w=%g\n",K[0],K[1],w); //bk
     2808//   if (w<0) {printf("negative w!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!\n");
     2809//             printf("K[0]=%g,K[1]=%g,w=%g\n",K[0],K[1],w);}
     2810
     2811
     2812}
     2813
    27892814/*}}}*/
    27902815/*FUNCTION Tria::Gradj {{{1*/
    … …  
    39153940                values[i]=solution[doflist[i]];
    39163941                if(isnan(values[i])) _error_("NaN found in solution vector");
    3917                 //if (values[i]<pow(10,-10))values[i]=pow(10,-10);
     3942        //      if (values[i]<pow(10,-10))values[i]=pow(10,-10); //correcting the water column to positive values
     3943 
    39183944        }
    39193945

  • issm/trunk/src/c/objects/IoModel.cpp

    r8538 r8549  
    211211        /*!Get thermal parameters: */
    212212        IoModelFetchData(&this->beta,iomodel_handle,"beta");
    213         IoModelFetchData(&this->hydro_gamma,iomodel_handle,"hydro_gamma");
    214         IoModelFetchData(&this->hydro_kn,iomodel_handle,"hydro_kn");
    215         IoModelFetchData(&this->hydro_p,iomodel_handle,"hydro_p");
    216         IoModelFetchData(&this->hydro_q,iomodel_handle,"hydro_q");
    217213        IoModelFetchData(&this->meltingpoint,iomodel_handle,"meltingpoint");
    218214        IoModelFetchData(&this->referencetemperature,iomodel_handle,"referencetemperature");
    … …  
    230226        IoModelFetchData(&this->numforcings,iomodel_handle,"numforcings");
    231227       
     228    /*!Get hydrology parameters: */
     229//      IoModelFetchData(&this->hydro_gamma,iomodel_handle,"hydro_gamma");
     230   IoModelFetchData(&this->hydro_kn,iomodel_handle,"hydro_kn");
     231   IoModelFetchData(&this->hydro_p,iomodel_handle,"hydro_p");
     232   IoModelFetchData(&this->hydro_q,iomodel_handle,"hydro_q");
     233        IoModelFetchData(&this->hydro_CR,iomodel_handle,"hydro_CR");
     234        IoModelFetchData(&this->hydro_n,iomodel_handle,"hydro_n");
    232235        /*qmu: */
    233236        if(this->qmu_analysis){

  • issm/trunk/src/c/objects/IoModel.h

    r8538 r8549  
    170170                /*thermal parameters: */
    171171                double beta;
    172                 double hydro_gamma;
    173                 double hydro_kn;
    174                 double hydro_p;
    175                 double hydro_q;
    176172                double meltingpoint;
    177173                double referencetemperature;
    … …  
    230226                int     loadcounter; //keep track of how many loads are being created in each analysis type
    231227                int     constraintcounter; //keep track of how many constraints are being created in each analysis type
    232                 /*}}}*/
     228
     229                 /*hydrology parameters: */
     230//               double hydro_gamma;
     231            double hydro_kn;
     232            double hydro_p;
     233                 double hydro_q;
     234       double hydro_CR;
     235                 double hydro_n;
     236                 /*}}}*/
    233237                /*Methods: {{{1*/
    234238                ~IoModel();

  • issm/trunk/src/c/objects/Materials/Matpar.cpp

    r8539 r8549  
    5050        matpar_mixed_layer_capacity=iomodel->mixed_layer_capacity;
    5151        matpar_thermal_exchange_velocity=iomodel->thermal_exchange_velocity;
    52         matpar_gamma=iomodel->hydro_gamma;
     52//      matpar_gamma=iomodel->hydro_gamma;
    5353        matpar_kn=iomodel->hydro_kn;
    5454
    … …  
    6969        this->hydro_p=iomodel->hydro_p;
    7070        this->hydro_q=iomodel->hydro_q;
    71 
     71   this->hydro_CR=iomodel->hydro_CR;
     72   this->hydro_n=iomodel->hydro_n;
    7273}
    7374/*}}}1*/
    … …  
    414415}
    415416/*}}}1*/
     417/*FUNCTION Matpar::GetHydroCR {{{1*/
     418double Matpar::GetHydroCR(){
     419           return hydro_CR;
     420}
     421/*}}}1*/
     422/*FUNCTION Matpar::GetHydroN {{{1*/
     423double Matpar::GetHydroN(){
     424           return hydro_n;
     425}
     426/*}}}1*/
    416427/*FUNCTION Matpar::TMeltingPoint {{{1*/
    417428double Matpar::TMeltingPoint(double pressure){

  • issm/trunk/src/c/objects/Materials/Matpar.h

    r8538 r8549  
    3333                double  hydro_p;
    3434                double  hydro_q;
    35 
     35      double  hydro_CR;
     36                double  hydro_n;
    3637        public:
    3738
    … …  
    8687                double GetHydroP();
    8788                double GetHydroQ();
     89      double GetHydroCR();
     90      double GetHydroN();
    8891                double TMeltingPoint(double pressure);
    8992                double PureIceEnthalpy(double pressure);

  • issm/trunk/src/c/solutions/hydrology_core.cpp

    r7643 r8549  
    4343        for(i=0;i<nsteps;i++){
    4444               
    45                 if(nsteps)_printf_(VerboseSolution(),"time step: %i/%i\n",i+1,nsteps);
     45                if(nsteps)_printf_(VerboseSolution(),"time step:%i/%i\n",i+1,nsteps);
    4646                time=(i+1)*dt;
    4747
     48//              printf("time step=%d=============================\n",i+1);//bk
    4849                /*call hydrology_core step: */
    4950                hydrology_core_step(femmodel,i,time);

  • issm/trunk/src/m/classes/version/7.6/model.m

    r8497 r8549  
    282282                 hydro_p=NaN;
    283283                 hydro_q=NaN;
    284                  hydro_gamma=0;
     284%                hydro_gamma=0;
    285285                 hydro_kn=0;
    286286                 
    … …  
    790790                         md.hydro_p=2;
    791791                         md.hydro_q=1;
    792                          md.hydro_gamma=1/(md.hydro_n*md.hydro_CR^md.hydro_p*(md.rho_water*md.g)^md.hydro_q);
     792%                        md.hydro_gamma=1/(md.hydro_n*md.hydro_CR^md.hydro_p*(md.rho_water*md.g)^md.hydro_q);
    793793                         md.hydro_kn=0;
    794794

  • issm/trunk/src/m/classes/version/7.7/model.m

    r8497 r8549  
    282282                 hydro_p=NaN;
    283283                 hydro_q=NaN;
    284                  hydro_gamma=0;
     284%                hydro_gamma=0;
    285285                 hydro_kn=0;
    286286                 
    … …  
    862862                         md.hydro_p=2;
    863863                         md.hydro_q=1;
    864                          md.hydro_gamma=1/(md.hydro_n*md.hydro_CR^md.hydro_p*(md.rho_water*md.g)^md.hydro_q);
     864%                        md.hydro_gamma=1/(md.hydro_n*md.hydro_CR^md.hydro_p*(md.rho_water*md.g)^md.hydro_q);
    865865                         md.hydro_kn=0;
    866866

  • issm/trunk/src/m/model/marshall.m

    r8541 r8549  
    200200%Thermal parameters
    201201WriteData(fid,md.beta,'Scalar','beta');
    202 WriteData(fid,md.hydro_gamma,'Scalar','hydro_gamma');
    203 WriteData(fid,md.hydro_kn,'Scalar','hydro_kn');
    204 WriteData(fid,md.hydro_p,'Scalar','hydro_p');
    205 WriteData(fid,md.hydro_q,'Scalar','hydro_q');
    206202WriteData(fid,md.meltingpoint,'Scalar','meltingpoint');
    207203WriteData(fid,md.referencetemperature,'Scalar','referencetemperature');
    … …  
    215211WriteData(fid,md.stabilize_constraints,'Integer','stabilize_constraints');
    216212
     213%Hydrology parameters
     214%WriteData(fid,md.hydro_gamma,'Scalar','hydro_gamma');
     215WriteData(fid,md.hydro_kn,'Scalar','hydro_kn');
     216WriteData(fid,md.hydro_p,'Scalar','hydro_p');
     217WriteData(fid,md.hydro_q,'Scalar','hydro_q');
     218WriteData(fid,md.hydro_CR,'Scalar','hydro_CR');
     219WriteData(fid,md.hydro_n,'Scalar','hydro_n');
     220
    217221%elements type
    218222WriteData(fid,md.ishutter,'Integer','ishutter');
Note: See TracChangeset for help on using the changeset viewer.