/////////////////////////////////////////////////////////////////////////////// // define, globales, fonctions perso (obligatoire et autres) ... // mohamed.berrada@upmc.fr /////////////////////////////////////////////////////////////////////////////// #include"../include/ncutil.h"//netcdf #define key_zco #include"../include/phy_cst_file.h" #include"../include/namelist.h" //#define PATH_NCFILES "/usr/home/mblod/neuro/nemo/nemo_opa/TRY/modipsl/config/GYRE/EXP00" #define PATH_NCFILES "../data_in" #include"../include/meshmask.h" // For PPCA /*double shfs_ta[NZ*NY*NX][NPCA]; double moy_ta[NZ*NY*NX]; double stdev_ta[NPCA]; #define shfs_ta( i, j, k ,n ) (shfs_ta[(i)*(NY*NZ)+(j)*(NZ)+(k)][(n)]) // shape functions #define moy_ta( i, j, k ) (moy_ta[(i)*(NY*NZ)+(j)*(NZ)+(k)]) // mean void load_shape_func (int argc, char *argv[]); // load PCA axes from file with mean void load_stdev_pac (int argc, char *argv[]); // load stdev from file void load_mean (int argc, char *argv[]); // load mean */ ////#define bmask( i, j ) (tmask( i, j ,0)) // ! elliptic equation is written at t-point int tniter[TA+TU]; // nbre d'iterations pour backward de solsor_dynspg_flt.h int ndiag_rst=0; // diagnostic restart if =1 //Arrays for neuler=1 (neuler=0::Euler 1st time step) #define tb_neuler1( i, j, k) (tb_neuler1[(k)*(NY*NX)+(j)*(NX)+(i)]) #define sb_neuler1( i, j, k) (sb_neuler1[(k)*(NY*NX)+(j)*(NX)+(i)]) #define ub_neuler1( i, j, k) (ub_neuler1[(k)*(NY*NX)+(j)*(NX)+(i)]) #define vb_neuler1( i, j, k) (vb_neuler1[(k)*(NY*NX)+(j)*(NX)+(i)]) #define sshb_neuler1( i, j) (sshb_neuler1[(j)*(NX)+(i)]) #define rotn_at_TU( i, j, k) (rotn_at_TU[(k)*(NY*NX)+(j)*(NX)+(i)]) #define hdivn_at_TU( i, j, k) (hdivn_at_TU[(k)*(NY*NX)+(j)*(NX)+(i)]) #define gcx_at_TU( i, j) (gcx_at_TU[(j)*(NX)+(i)]) double tb_neuler1[NZ*NY*NX]; // double sb_neuler1[NZ*NY*NX]; // double ub_neuler1[NZ*NY*NX]; // double vb_neuler1[NZ*NY*NX]; // double sshb_neuler1[NY*NX]; // double rotn_at_TU[NZ*NY*NX]; // double hdivn_at_TU[NZ*NY*NX]; // double gcx_at_TU[NY*NX]; // void true_target_in_tab(int argc, char *argv[]); // charge YS_*a_c(0,i,j,k) dans true_*a_c(i,j,k) void normsup_alla_c(); // diagnostic of minimization double gdsr[NZ]; int nksr; void xtraqsr_init();// init of solar radiation penetration # define rdttra( k ) rdt double r2dt; double atfp1 = 1. - 2. * atfp;// !: asselin time filter coeff. (atfp1= 1-2*atfp) //------->zdfmxl double avt_c = 5.e-4; // ! Kz criterion for the turbocline depth double rho_c = 0.01;// ! density criterion for mixed layer depth double ppgphi0 = 29.; // latitude of first raw column T-point //latitude for the Coriolis or Beta parameter //mode de sauvegarde (==0 3D, ==1 2D surface) short savemode=0; int n_zdfexp=3; int jphgr_msh=5;//pour ff :beta-plane with regular grid-spacing and rotated domain (GYRE config) // a faire /*nyear = ndastp / 10000 : current year nmonth : current month of the year nyear nday_year : current day of the year nyea ndastp : = nyear*10000 + nmonth*100 + nday*/ void cal_ff(); void xdisplay() ; void xistate_gyre(); void xistate_yao_file(char * pth_file_yao); // a faire void xdom_init(); void xflt_rst(); void xsolmat_init(); void xistate_init(int argc, char *argv[]); void xrst_save(int argc, char *argv[]); ////////////////////////////////////////////////////////////// // Les fonctions OBLIGATOIRES //_____________________________________________________________________________ void appli_start (int argc, char *argv[]) { // permet si besoin de prendre la main des le debut de l'appl; printf("////////////////////////////////////////////////////////////////////////\n"); printf("// NEMO/YAO PROJECT //\n"); printf("// M. Berrada 02-2009 //\n"); printf("// LOCEAN-IPSL.UPMC (Paris 6) //\n"); printf("//====================================================================//\n"); printf("\n"); #ifdef K_SOLSORYAO printf("option SOLSORYAO enable\n"); #else printf("option SOLSORYAO disable\n"); #endif #ifdef K_OPTIMORDER printf("option OPTIMORDER enable\n"); #else printf("option OPTIMORDER disable\n"); #endif //Vérification du type de réel utilisé //A modifier - pas génial #ifdef YFLOAT #ifndef YYFLOAT printf("Incoherent real type between O_REAL option and ncutil.h\n"); exit(1); #endif #endif #ifdef YDOUBLE #ifndef YYDOUBLE printf("Incoherent real type between O_REAL option and ncutil.h\n"); exit(1); #endif #endif for (int i=0;i appel via loadstate // =YIO_SAVESTATE => appel via savestate // =YIO_LOADOBS => appel via loadobs // =YIO_OUTOOBS => appel via outoobs if(indic==YIO_SAVESTATE && !strcmp(nmmod,"tb") && savemode){ if(kaxe==0 ) return(1); else return(0); } /* if(indic==YIO_OUTOOBS){ if(iaxe==9 && jaxe==9 ) return(1); else return(0); }*/ return(1); } //======================== Fin des foncts obligatoires ========================= //============================================================================== //Les autres fonctions //_____________________________________________________________________________ void xdisplay() { normsup_alla_c(); } //_____________________________________________________________________________ #define hu( i, j ) (hu[(j)*(NX)+(i)]) #define hv( i, j ) (hv[(j)*(NX)+(i)]) #define hur( i, j ) (hur[(j)*(NX)+(i)]) #define hvr( i, j ) (hvr[(j)*(NX)+(i)]) double hu[NY*NX]; // xdom_init double hv[NY*NX]; // xdom_init double hur[NY*NX]; // xdom_init double hvr[NY*NX]; // xdom_init void xdom_init(){ for(int j=0;jNZ-1) nksr=NZ-1; } } //_____________________________________________________________________________ void cal_ff() //domhgr.F90 Coriolis parameter { if(jphgr_msh==5){// beta-plane and rotated domain (gyre configuration) double zbeta = 2.*omega*cos(rad*ppgphi0)/ra ;// beta at latitude ppgphi0 double zphi0 = 15.;//latitude of the first row F-points double zf0 = 2.*omega*sin(rad*zphi0);// compute f0 1st point south for(int i=0;i