New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 4171 – NEMO

Changeset 4171


Ignore:
Timestamp:
2013-11-08T18:29:37+01:00 (10 years ago)
Author:
rfurner
Message:

merging changes from 3987:4141 of dev_r3987_UKMO4_OBS

Location:
branches/2013/dev_UKMO_2013
Files:
1 deleted
26 edited
6 copied

Legend:

Unmodified
Added
Removed

  • branches/2013/dev_UKMO_2013/DOC/TexFiles/Chapters/Chap_DIA.tex

    r3976 r4171  
    362362\subsubsection{Use of Groups} 
    363363 
    364 Groups can be used for 2 purposes. Firstly, the group can be used to define common attributes to be shared by the elements of the group through the inheritance. In the following example, we define a group of field that will share a common grid ''grid\_T\_2D''. Note that for the field ''toce'', we overwrite the grid definition inherited from the group by ''grid\_T\_3D''. 
     364Groups can be used for 2 purposes. Firstly, the group can be used to define common attributes to be shared by the elements of the group through inheritance. In the following example, we define a group of field that will share a common grid ''grid\_T\_2D''. Note that for the field ''toce'', we overwrite the grid definition inherited from the group by ''grid\_T\_3D''. 
    365365\vspace{-20pt} 
    366366\begin{alltt}  {{\scriptsize 
    … …  
    386386\end{verbatim} 
    387387}}\end{alltt}  
    388 that can be directly include in a file through the following syntax: 
     388that can be directly included in a file through the following syntax: 
    389389\vspace{-20pt} 
    390390\begin{alltt}  {{\scriptsize 
    … …  
    466466\end{verbatim} 
    467467}}\end{alltt}  
    468 However it is often very convienent to define the file name with the name of the experience, the output file frequency and the date of the beginning and the end of the simulation (which are informations stored either in the namelist or in the XML file). To do so, we added the following rule: if the id of the tag file is ''fileN''(where N = 1 to 99) or one of the predefined section or mooring (see next subsection), the following part of the name and the name\_suffix (that can be inherited) will be automatically replaced by:\\ 
     468However it is often very convienent to define the file name with the name of the experiment, the output file frequency and the date of the beginning and the end of the simulation (which are informations stored either in the namelist or in the XML file). To do so, we added the following rule: if the id of the tag file is ''fileN''(where N = 1 to 99) or one of the predefined sections or moorings (see next subsection), the following part of the name and the name\_suffix (that can be inherited) will be automatically replaced by:\\ 
    469469\\ 
    470470\begin{tabular}{|p{4cm}|p{8cm}|} 
    … …  
    474474   \hline 
    475475   \centering @expname@ & 
    476    the experience name (from cn\_exp in the namelist) \\ 
     476   the experiment name (from cn\_exp in the namelist) \\ 
    477477   \hline 
    478478   \centering @freq@ & 
    … …  
    590590   file\_definition &  
    591591   encapsulates the definition of all the files that will be outputted & 
    592    enabled, min\_digits, name, name\_suffix, output\_level, split\_format, split\_freq, sync\_freq, type, src & 
     592   enabled, min\_digits, name, name\_suffix, output\_level, split\_freq\_format, split\_freq, sync\_freq, type, src & 
    593593   context &  
    594594   file or file\_group \\ 
    … …  
    596596   file\_group &  
    597597   encapsulates a group of files that will be outputted & 
    598    enabled, description, id, min\_digits, name, name\_suffix, output\_freq, output\_level, split\_format, split\_freq, sync\_freq, type, src & 
     598   enabled, description, id, min\_digits, name, name\_suffix, output\_freq, output\_level, split\_freq\_format, split\_freq, sync\_freq, type, src & 
    599599   file\_definition, file\_group &  
    600600   file or file\_group \\ 
    … …  
    602602   file &  
    603603   define the contents of a file to be outputted & 
    604    enabled, description, id, min\_digits, name, name\_suffix, output\_freq, output\_level, split\_format, split\_freq, sync\_freq, type, src & 
     604   enabled, description, id, min\_digits, name, name\_suffix, output\_freq, output\_level, split\_freq\_format, split\_freq, sync\_freq, type, src & 
    605605   file\_definition, file\_group &  
    606606   field \\ 
    … …  
    775775   field family \\  
    776776   \hline    
    777    split\_format &  
    778    date format used in the name of splitted output files. can be spécified using the following syntaxe: \%y, \%mo, \%d, \%h \%mi and \%s &  
    779    split\_format= "\%yy\%mom\%dd" &  
     777   split\_freq &  
     778   frequency at which to temporally split output files. Units can be ts (timestep), y, mo, d, h, mi, s. Useful for long runs to prevent over-sized output files.&  
     779   split\_freq="1mo" &  
    780780   file family \\  
    781781   \hline    
    782    split\_freq &  
    783    split output files frequency. units can be ts (timestep), y, mo, d, h, mi, s. &  
    784    split\_freq="1mo" &  
     782   split\_freq\-\_format &  
     783   date format used in the name of temporally split output files. Can be specified  
     784   using the following syntaxes: \%y, \%mo, \%d, \%h \%mi and \%s &  
     785   split\_freq\_format= "\%y\%mo\%d" &  
    785786   file family \\  
    786787   \hline    
    … …  
    812813   \hline    
    813814   type (1)&  
    814    specify if the output files must be split (multiple\_file) or not (one\_file) &  
     815   specify if the output files are to be split spatially (multiple\_file) or not (one\_file) &  
    815816   type="multiple\_file" &  
    816817   file familly \\  

  • branches/2013/dev_UKMO_2013/DOC/TexFiles/Chapters/Chap_OBS.tex

    r3294 r4171  
    55\label{OBS} 
    66 
    7 Authors: D. Lea, M. Martin, K. Mogensen, A. Vidard, A. Weaver...   % do we keep that ? 
     7Authors: D. Lea, M. Martin, K. Mogensen, A. Vidard, A. Weaver, A. Ryan, ...   % do we keep that ? 
    88 
    99\minitoc 
    … …  
    4242where to obtain data and how to setup the namelist. Section~\ref{OBS_details} introduces some 
    4343more technical details of the different observation types used and also shows a more complete 
    44 namelist. Section~\ref{OBS_theory} introduces some of the theoretical aspects of the 
    45 observation operator including interpolation methods and running on multiple processors. 
    46 Section~\ref{OBS_obsutils} introduces some utilities to help working with the files produced 
    47 by the OBS code. 
     44namelist. Section~\ref{OBS_theory} introduces some of the theoretical aspects of the observation 
     45operator including interpolation methods and running on multiple processors. 
     46Section~\ref{OBS_ooo} describes the offline observation operator code. 
     47Section~\ref{OBS_obsutils} introduces some utilities to help working with the files 
     48produced by the OBS code. 
    4849 
    4950% ================================================================ 
    … …  
    786787\newpage 
    787788 
     789% ================================================================ 
     790% Offline observation operator documentation 
     791% ================================================================ 
     792 
     793%\usepackage{framed} 
     794 
     795\section{Offline observation operator} 
     796\label{OBS_ooo} 
     797 
     798\subsection{Concept} 
     799 
     800The obs oper maps model variables to observation space. It is possible to apply this mapping 
     801without running the model. The software which performs this functionality is known as the 
     802\textbf{offline obs oper}. The obs oper is divided into three stages. An initialisation phase, 
     803an interpolation phase and an output phase. The implementation of which is outlined in the 
     804previous sections. During the interpolation phase the offline obs oper populates the model 
     805arrays by reading saved model fields from disk. 
     806 
     807There are two ways of exploiting this offline capacity. The first is to mimic the behaviour of 
     808the online system by supplying model fields at regular intervals between the start and the end 
     809of the run. This approach results in a single model counterpart per observation. This kind of 
     810usage produces feedback files the same file format as the online obs oper.  
     811The second is to take advantage of the offline setting in which multiple model counterparts can 
     812be calculated per observation. In this case it is possible to consider all forecasts verifying 
     813at the same time. By forecast, I mean any method which produces an estimate of physical reality 
     814which is not an observed value. In the case of class 4 files this means forecasts, analyses, persisted 
     815analyses and climatological values verifying at the same time. Although the class 4 file format 
     816doesn't account for multiple ensemble members or multiple experiments per observation, it is possible 
     817to include these components in the same or multiple files. 
     818 
     819%-------------------------------------------------------------------------------------------------------- 
     820% offline_oper.exe  
     821%-------------------------------------------------------------------------------------------------------- 
     822 
     823\subsection{Using the offline observation operator} 
     824 
     825\subsubsection{Building} 
     826 
     827In addition to \emph{OPA\_SRC} the offline obs oper requires the inclusion 
     828of the \emph{OOO\_SRC} directory. \emph{OOO\_SRC} contains a replacement \textbf{nemo.f90} and 
     829\textbf{nemogcm.F90} which overwrites the resultant \textbf{nemo.exe}. This is the approach taken 
     830by \emph{SAS\_SRC} and \emph{OFF\_SRC}. 
     831 
     832%-------------------------------------------------------------------------------------------------------- 
     833% Running  
     834%-------------------------------------------------------------------------------------------------------- 
     835\subsubsection{Running} 
     836 
     837The simplest way to use the executable is to edit and append the \textbf{ooo.nml} namelist to 
     838a full NEMO namelist and then to run the executable as if it were nemo.exe.  
     839 
     840\subsubsection{Quick script} 
     841 
     842A useful Python utility to control the namelist options can be found in \textbf{OBSTOOLS/OOO}. The 
     843functions which locate model fields and observation files can be manually specified. The package 
     844can be installed by appropriate use of the included setup.py script. 
     845 
     846Documentation can be auto-generated by Sphinx by running \emph{make html} in the \textbf{doc} directory. 
     847 
     848%-------------------------------------------------------------------------------------------------------- 
     849% Configuration section 
     850%-------------------------------------------------------------------------------------------------------- 
     851\subsection{Configuring the offline observation operator} 
     852The observation files and settings understood by \textbf{namobs} have been outlined in the online 
     853obs oper section. In addition there are two further namelists wich control the operation of the offline 
     854obs oper. \textbf{namooo} which controls the input model fields and \textbf{namcl4} which controls the 
     855production of class 4 files.  
     856 
     857\subsubsection{Single field} 
     858 
     859In offline mode model arrays are populated at appropriate time steps via input files. 
     860At present, \textbf{tsn} and \textbf{sshn} are populated by the default read routines.  
     861These routines will be expanded upon in future versions to allow the specification of any 
     862model variable. As such, input files must be global versions of the model domain with 
     863\textbf{votemper}, \textbf{vosaline} and optionally \textbf{sshn} present. 
     864 
     865For each field read there must be an entry in the \textbf{namooo} namelist specifying the 
     866name of the file to read and the index along the \emph{time\_counter}. For example, to 
     867read the second time counter from a single file the namelist would be. 
     868 
     869\begin{alltt} 
     870\tiny 
     871\begin{verbatim}  
     872!---------------------------------------------------------------------- 
     873!       namooo Offline obs_oper namelist 
     874!---------------------------------------------------------------------- 
     875!   ooo_files    specifies the files containing the model counterpart 
     876!   nn_ooo_idx   specifies the time_counter index within the model file 
     877&namooo 
     878   ooo_files = "foo.nc" 
     879   nn_ooo_idx = 2 
     880/ 
     881\end{verbatim}  
     882\end{alltt} 
     883 
     884\subsubsection{Multiple fields per run} 
     885 
     886Model field iteration is controlled via \textbf{nn\_ooo\_freq} which specifies 
     887the number of model steps at which the next field gets read. For example, if 
     88812 hourly fields are to be interpolated in a setup where 288 steps equals 24 hours. 
     889 
     890\begin{alltt} 
     891\tiny 
     892\begin{verbatim}  
     893!---------------------------------------------------------------------- 
     894!       namooo Offline obs_oper namelist 
     895!---------------------------------------------------------------------- 
     896!   ooo_files    specifies the files containing the model counterpart 
     897!   nn_ooo_idx   specifies the time_counter index within the model file 
     898!   nn_ooo_freq  specifies number of time steps between read operations 
     899&namooo 
     900   ooo_files = "foo.nc" "foo.nc" 
     901   nn_ooo_idx = 1 2 
     902   nn_ooo_freq = 144 
     903/ 
     904\end{verbatim}  
     905\end{alltt} 
     906 
     907The above namelist will result in feedback files whose first 12 hours contain 
     908the first field of foo.nc and the second 12 hours contain the second field. 
     909 
     910%\begin{framed} 
     911\textbf{Note} Missing files can be denoted as "nofile". 
     912%\end{framed} 
     913 
     914It is easy to see how a collection of fields taken fron a number of files  
     915at different indices can be combined at a particular frequency in time to  
     916generate a pseudo model evolution. As long as all that is needed is a single 
     917model counterpart at a regular interval then namooo is all that needs to 
     918be edited. However, a far more interesting approach can be taken in which 
     919multiple forecasts, analyses, persisted analyses and climatologies are 
     920considered against the same set of observations. For this a slightly more 
     921complicated approach is needed. It is referred to as \emph{Class 4} since 
     922it is the fourth metric defined by the GODAE intercomparison project. 
     923 
     924%-------------------------------------------------------------------------------------------------------- 
     925% Class 4 file section 
     926%-------------------------------------------------------------------------------------------------------- 
     927\subsubsection{Multiple model counterparts per observation a.k.a Class 4} 
     928 
     929A generalisation of feedback files to allow multiple model components per observation. For a single 
     930observation, as well as previous forecasts verifying at the same time there are also analyses, persisted 
     931analyses and climatologies.  
     932 
     933 
     934The above namelist performs two basic functions. It organises the fields 
     935given in \textbf{namooo} into groups so that observations can be matched 
     936up multiple times. It also controls the metadata and the output variable 
     937of the class 4 file when a write routine is called. 
     938 
     939%\begin{framed} 
     940\textbf{Note: ln\_cl4} must be set to \emph{.TRUE.} in \textbf{namobs}  
     941to use class 4 outputs. 
     942%\end{framed} 
     943 
     944\subsubsection{Class 4 naming convention} 
     945 
     946The standard class 4 file naming convention is as follows. 
     947 
     948\noindent 
     949\linebreak 
     950\textbf{\$\{prefix\}\_\$\{yyyymmdd\}\_\$\{sys\}\_\$\{cfg\}\_\$\{vn\}\_\$\{kind\}\_\$\{nproc\}.nc} 
     951 
     952\noindent 
     953\linebreak 
     954Much of the namelist is devoted to specifying this convention. The 
     955following namelist settings control the elements of the output 
     956file names. Each should be specified as a single string of character data. 
     957 
     958\begin{description} 
     959\item[cl4\_prefix] 
     960Prefix for class 4 files e.g. class4 
     961\item[cl4\_date] 
     962YYYYMMDD validity date 
     963\item[cl4\_sys] 
     964The name of the class 4 model system e.g. FOAM 
     965\item[cl4\_cfg] 
     966The name of the class 4 model configuration e.g. orca025 
     967\item[cl4\_vn] 
     968The name of the class 4 model version e.g. 12.0 
     969\end{description} 
     970 
     971\noindent 
     972The kind is specified by the observation type internally to the obs oper. The processor 
     973number is specified internally in NEMO.  
     974 
     975\subsubsection{Class 4 file global attributes} 
     976 
     977Global attributes necessary to fulfill the class 4 file definition. These 
     978are also useful pieces of information when collaborating with external 
     979partners. 
     980 
     981\begin{description} 
     982\item[cl4\_contact] 
     983Contact email for class 4 files. 
     984\item[cl4\_inst] 
     985The name of the producers institution. 
     986\item[cl4\_cfg] 
     987The name of the class 4 model configuration e.g. orca025 
     988\item[cl4\_vn] 
     989The name of the class 4 model version e.g. 12.0 
     990\end{description} 
     991 
     992\noindent 
     993The obs\_type, 
     994creation date and validity time are specified internally to the obs oper. 
     995 
     996\subsubsection{Class 4 model counterpart configuration} 
     997 
     998As seen previously it is possible to perform a single sweep of the 
     999obs oper and specify a collection of model fields equally spaced  
     1000along that sweep. In the class 4 case the single sweep is replaced 
     1001with multiple sweeps and a certain ammount of book keeping is 
     1002needed to ensure each model counterpart makes its way to the  
     1003correct piece of memory in the output files. 
     1004 
     1005\noindent 
     1006\linebreak 
     1007In terms of book keeping, the offline obs oper needs to know how many 
     1008full sweeps need to be performed. This is specified via the  
     1009\textbf{cl4\_match\_len} variable and is the total number of model 
     1010counterparts per observation. For example, a 3 forecasts plus 3 persistence 
     1011fields plus an analysis field would be 7 counterparts per observation. 
     1012 
     1013\begin{alltt} 
     1014\tiny 
     1015\begin{verbatim} 
     1016   cl4_match_len = 7 
     1017\end{verbatim} 
     1018\end{alltt} 
     1019 
     1020Then to correctly allocate a class 4 file the forecast axis must be defined. This 
     1021is controlled via \textbf{cl4\_fcst\_len}, which in out above example would be 3. 
     1022 
     1023\begin{alltt} 
     1024\tiny 
     1025\begin{verbatim} 
     1026   cl4_fcst_len = 3 
     1027\end{verbatim} 
     1028\end{alltt} 
     1029 
     1030Then for each model field it is necessary to designate what class 4 variable and 
     1031index along the forecast dimension the model counterpart should be stored in the 
     1032output file. As well as a value for that lead time in hours, this will be useful 
     1033when interpreting the data afterwards.  
     1034 
     1035\begin{alltt} 
     1036\tiny 
     1037\begin{verbatim} 
     1038   cl4_vars = "forecast" "forecast" "forecast" "persistence" "persistence" 
     1039              "persistence" "best_estimate" 
     1040   cl4_fcst_idx = 1 2 3 1 2 3 1 
     1041   cl4_leadtime = 12 36 60  
     1042\end{verbatim} 
     1043\end{alltt} 
     1044 
     1045In terms of files and indices of fields inside each file the class 4 approach 
     1046makes use of the \textbf{namooo} namelist. If our fields are in separate files 
     1047with a single field per file our example inputs will be specified. 
     1048 
     1049\begin{alltt} 
     1050\tiny 
     1051\begin{verbatim} 
     1052   ooo_files = "F.1.nc" "F.2.nc" "F.3.nc" "P.1.nc" "P.2.nc" "P.3.nc" "A.1.nc" 
     1053   nn_ooo_idx = 1 1 1 1 1 1 1 
     1054\end{verbatim} 
     1055\end{alltt} 
     1056 
     1057When we combine all of the naming conventions, global attributes and i/o instructions 
     1058the class 4 namelist becomes. 
     1059 
     1060\begin{alltt} 
     1061\tiny 
     1062\begin{verbatim} 
     1063!---------------------------------------------------------------------- 
     1064!       namooo Offline obs_oper namelist 
     1065!---------------------------------------------------------------------- 
     1066!   ooo_files    specifies the files containing the model counterpart 
     1067!   nn_ooo_idx   specifies the time_counter index within the model file 
     1068!   nn_ooo_freq  specifies number of time steps between read operations 
     1069&namooo 
     1070   ooo_files = "F.1.nc" "F.2.nc" "F.3.nc" "P.1.nc" "P.2.nc" "P.3.nc" "A.1.nc" 
     1071   nn_ooo_idx = 1 1 1 1 1 1 1 
     1072/ 
     1073!---------------------------------------------------------------------- 
     1074!       namcl4 Offline obs_oper class 4 namelist 
     1075!---------------------------------------------------------------------- 
     1076! 
     1077!  Naming convention 
     1078!  ----------------- 
     1079!  cl4_prefix    specifies the output file prefix 
     1080!  cl4_date      specifies the output file validity date 
     1081!  cl4_sys       specifies the model counterpart system 
     1082!  cl4_cfg       specifies the model counterpart configuration 
     1083!  cl4_vn        specifies the model counterpart version 
     1084!  cl4_inst      specifies the model counterpart institute 
     1085!  cl4_contact   specifies the file producers contact details 
     1086! 
     1087!  I/O specification 
     1088!  ----------------- 
     1089!  cl4_vars      specifies the names of the output file netcdf variable 
     1090!  cl4_fcst_idx  specifies output file forecast index 
     1091!  cl4_fcst_len  specifies forecast axis length 
     1092!  cl4_match_len specifies number of unique matches per observation 
     1093!  cl4_leadtime  specifies the forecast axis lead time  
     1094! 
     1095&namcl4 
     1096   cl4_match_len = 7 
     1097   cl4_fcst_len = 3 
     1098   cl4_fcst_idx = 1 2 3 1 2 3 1 
     1099   cl4_vars = "forecast" "forecast" "forecast" "persistence" "persistence" 
     1100              "persistence" "best_estimate" 
     1101   cl4_leadtime = 12 36 60 
     1102   cl4_prefix = "class4" 
     1103   cl4_date = "20130101" 
     1104   cl4_vn = "12.0" 
     1105   cl4_sys = "FOAM" 
     1106   cl4_cfg = "AMM7" 
     1107   cl4_contact = "example@example.com" 
     1108   cl4_inst = "UK Met Office" 
     1109/ 
     1110\end{verbatim} 
     1111\end{alltt} 
     1112 
     1113\subsubsection{Climatology interpolation} 
     1114 
     1115The climatological counterpart is generated at the start of the run by restarting  
     1116the model from climatology through appropriate use of \textbf{namtsd}. To override 
     1117the offline observation operator read routine and to take advantage of the restart 
     1118settings, specify the first entry in \textbf{cl4\_vars} as "climatology". This will then 
     1119pipe the restart from climatology into the output class 4 file. As in every other 
     1120class 4 matchup the input file, input index and output index must be specified. 
     1121These can be replaced with dummy data since they are not used but they must be 
     1122present to cycle through the matchups correctly.  
     1123 
     1124\subsection{Advanced usage} 
     1125 
     1126In certain cases it may be desirable to include both multiple model fields per 
     1127observation window with multiple match ups per observation. This can be achieved 
     1128by specifying \textbf{nn\_ooo\_freq} as well as the class 4 settings. Care must 
     1129be taken in generating the ooo\_files list such that the files are arranged into 
     1130consecutive blocks of single match ups. For example, 2 forecast fields  
     1131of 12 hourly data would result in 4 separate read operations but only 2 write 
     1132operations, 1 per forecast. 
     1133 
     1134\begin{alltt} 
     1135\tiny 
     1136\begin{verbatim} 
     1137   ooo_files = "F1.nc" "F1.nc" "F2.nc" "F2.nc" 
     1138... 
     1139   cl4_fcst_idx = 1 2 
     1140\end{verbatim} 
     1141\end{alltt} 
     1142 
     1143The above notation reveals the internal split between match up iterators and file 
     1144iterators. This technique has not been used before so experimentation is needed 
     1145before results can be trusted. 
     1146 
     1147 
     1148 
     1149 
     1150\newpage 
     1151 
    7881152\section{Observation Utilities} 
    7891153\label{OBS_obsutils} 
    … …  
    8011165handling observation files and the feedback file output from the NEMO observation operator. 
    8021166The utilities are as follows 
     1167 
     1168\subsubsection{c4comb} 
     1169 
     1170The program c4comb combines multiple class 4 files produced by individual processors in an 
     1171MPI run of NEMO offline obs\_oper into a single class 4 file. The program is called in the following way: 
     1172 
     1173\begin{alltt} 
     1174\footnotesize 
     1175\begin{verbatim} 
     1176c4comb.exe outputfile inputfile1 inputfile2 ... 
     1177\end{verbatim} 
     1178\end{alltt} 
    8031179 
    8041180\subsubsection{corio2fb} 

  • branches/2013/dev_UKMO_2013/NEMOGCM/NEMO/NST_SRC/agrif_opa_sponge.F90

    r3918 r4171  
    185185      INTEGER  :: ji,jj,jk 
    186186      INTEGER  :: ispongearea, ilci, ilcj 
    187       REAL(wp) :: z1spongearea 
    188       REAL(wp), POINTER, DIMENSION(:,:) :: zlocalviscsponge 
     187      LOGICAL  :: ll_spdone 
     188      REAL(wp) :: z1spongearea, zramp 
     189      REAL(wp), POINTER, DIMENSION(:,:) :: ztabramp 
    189190 
    190191#if defined SPONGE || defined SPONGE_TOP 
    191  
    192       CALL wrk_alloc( jpi, jpj, zlocalviscsponge ) 
    193  
    194       ispongearea  = 2 + 2 * Agrif_irhox() 
    195       ilci = nlci - ispongearea 
    196       ilcj = nlcj - ispongearea  
    197       z1spongearea = 1._wp / REAL( ispongearea - 2 ) 
    198       spbtr2(:,:) = 1. / ( e1t(:,:) * e2t(:,:) ) 
     192      ll_spdone=.TRUE. 
     193      IF (( .NOT. spongedoneT ).OR.( .NOT. spongedoneU )) THEN 
     194         ! Define ramp from boundaries towards domain interior 
     195         ! at T-points 
     196         ! Store it in ztabramp 
     197         ll_spdone=.FALSE. 
     198 
     199         CALL wrk_alloc( jpi, jpj, ztabramp ) 
     200 
     201         ispongearea  = 2 + 2 * Agrif_irhox() 
     202         ilci = nlci - ispongearea 
     203         ilcj = nlcj - ispongearea  
     204         z1spongearea = 1._wp / REAL( ispongearea - 2 ) 
     205         spbtr2(:,:) = 1. / ( e1t(:,:) * e2t(:,:) ) 
     206 
     207         ztabramp(:,:) = 0. 
     208 
     209         IF( (nbondi == -1) .OR. (nbondi == 2) ) THEN 
     210            DO jj = 1, jpj 
     211               IF ( umask(2,jj,1) == 1._wp ) THEN 
     212                 DO ji = 2, ispongearea                   
     213                    ztabramp(ji,jj) = ( ispongearea-ji ) * z1spongearea 
     214                 END DO 
     215               ENDIF 
     216            ENDDO 
     217         ENDIF 
     218 
     219         IF( (nbondi == 1) .OR. (nbondi == 2) ) THEN 
     220            DO jj = 1, jpj 
     221               IF ( umask(nlci-2,jj,1) == 1._wp ) THEN 
     222                  DO ji = ilci+1,nlci-1 
     223                     zramp = (ji - (ilci+1) ) * z1spongearea 
     224                     ztabramp(ji,jj) = MAX( ztabramp(ji,jj), zramp ) 
     225                  ENDDO 
     226               ENDIF 
     227            ENDDO 
     228         ENDIF 
     229 
     230         IF( (nbondj == -1) .OR. (nbondj == 2) ) THEN 
     231            DO ji = 1, jpi 
     232               IF ( vmask(ji,2,1) == 1._wp ) THEN 
     233                  DO jj = 2, ispongearea 
     234                     zramp = ( ispongearea-jj ) * z1spongearea 
     235                     ztabramp(ji,jj) = MAX( ztabramp(ji,jj), zramp ) 
     236                  END DO 
     237               ENDIF 
     238            ENDDO 
     239         ENDIF 
     240 
     241         IF( (nbondj == 1) .OR. (nbondj == 2) ) THEN 
     242            DO ji = 1, jpi 
     243               IF ( vmask(ji,nlcj-2,1) == 1._wp ) THEN 
     244                  DO jj = ilcj+1,nlcj-1 
     245                     zramp = (jj - (ilcj+1) ) * z1spongearea 
     246                     ztabramp(ji,jj) = MAX( ztabramp(ji,jj), zramp ) 
     247                  END DO 
     248               ENDIF 
     249            ENDDO 
     250         ENDIF 
     251 
     252      ENDIF 
    199253 
    200254      ! Tracers 
    201255      IF( .NOT. spongedoneT ) THEN 
    202          zlocalviscsponge(:,:) = 0. 
    203256         spe1ur(:,:) = 0. 
    204257         spe2vr(:,:) = 0. 
    205258 
    206259         IF( (nbondi == -1) .OR. (nbondi == 2) ) THEN 
    207             DO ji = 2, ispongearea 
    208                zlocalviscsponge(ji,:) = visc_tra * ( ispongearea-ji ) * z1spongearea 
    209             ENDDO 
    210             spe1ur(2:ispongearea-1,:      ) = 0.5 * ( zlocalviscsponge(2:ispongearea-1,:      )   & 
    211                &                         +            zlocalviscsponge(3:ispongearea  ,:      ) ) & 
    212                &                         * e2u(2:ispongearea-1,:      ) / e1u(2:ispongearea-1,:      ) 
    213             spe2vr(2:ispongearea  ,1:jpjm1) = 0.5 * ( zlocalviscsponge(2:ispongearea  ,1:jpjm1)   & 
    214                &                         +            zlocalviscsponge(2:ispongearea,2  :jpj  ) ) & 
    215                &                         * e1v(2:ispongearea  ,1:jpjm1) / e2v(2:ispongearea  ,1:jpjm1) 
     260            spe1ur(2:ispongearea-1,:       ) = visc_tra                                        & 
     261               &                             *    0.5 * (  ztabramp(2:ispongearea-1,:      )   & 
     262               &                                         + ztabramp(3:ispongearea  ,:      ) ) & 
     263               &                             * e2u(2:ispongearea-1,:) / e1u(2:ispongearea-1,:) 
     264 
     265            spe2vr(2:ispongearea  ,1:jpjm1 ) = visc_tra                                        & 
     266               &                             *    0.5 * (  ztabramp(2:ispongearea  ,1:jpjm1)   & 
     267               &                                         + ztabramp(2:ispongearea,2  :jpj  ) ) & 
     268               &                             * e1v(2:ispongearea,1:jpjm1) / e2v(2:ispongearea,1:jpjm1) 
    216269         ENDIF 
    217270 
    218271         IF( (nbondi == 1) .OR. (nbondi == 2) ) THEN 
    219             DO ji = ilci+1,nlci-1 
    220                zlocalviscsponge(ji,:) = visc_tra * (ji - (ilci+1) ) * z1spongearea 
    221             ENDDO 
    222    
    223             spe1ur(ilci+1:nlci-2,:      ) = 0.5 * (  zlocalviscsponge(ilci+1:nlci-2,:)    &  
    224                &                          +          zlocalviscsponge(ilci+2:nlci-1,:) )  & 
    225                &                          * e2u(ilci+1:nlci-2,:) / e1u(ilci+1:nlci-2,:) 
    226  
    227             spe2vr(ilci+1:nlci-1,1:jpjm1) = 0.5 * (  zlocalviscsponge(ilci+1:nlci-1,1:jpjm1)    &  
    228                &                            +        zlocalviscsponge(ilci+1:nlci-1,2:jpj  )  ) &  
    229                &                                   * e1v(ilci+1:nlci-1,1:jpjm1) / e2v(ilci+1:nlci-1,1:jpjm1) 
     272            spe1ur(ilci+1:nlci-2,:        ) = visc_tra                                   & 
     273               &                            * 0.5 * (  ztabramp(ilci+1:nlci-2,:      )   &  
     274               &                                     + ztabramp(ilci+2:nlci-1,:      ) ) & 
     275               &                            * e2u(ilci+1:nlci-2,:) / e1u(ilci+1:nlci-2,:) 
     276 
     277            spe2vr(ilci+1:nlci-1,1:jpjm1  )  = visc_tra                                  & 
     278               &                            * 0.5 * (  ztabramp(ilci+1:nlci-1,1:jpjm1)   &  
     279               &                                     + ztabramp(ilci+1:nlci-1,2:jpj  ) ) &  
     280               &                            * e1v(ilci+1:nlci-1,1:jpjm1) / e2v(ilci+1:nlci-1,1:jpjm1) 
    230281         ENDIF 
    231282 
    232283         IF( (nbondj == -1) .OR. (nbondj == 2) ) THEN 
    233             DO jj = 2, ispongearea 
    234                zlocalviscsponge(:,jj) = visc_tra * ( ispongearea-jj ) * z1spongearea 
    235             ENDDO 
    236             spe1ur(1:jpim1,2:ispongearea  ) = 0.5 * ( zlocalviscsponge(1:jpim1,2:ispongearea  ) &  
    237                &                            +         zlocalviscsponge(2:jpi  ,2:ispongearea) ) & 
     284            spe1ur(1:jpim1,2:ispongearea  ) = visc_tra                                     & 
     285               &                            * 0.5 * (  ztabramp(1:jpim1,2:ispongearea  )   &  
     286               &                                     + ztabramp(2:jpi  ,2:ispongearea  ) ) & 
    238287               &                            * e2u(1:jpim1,2:ispongearea) / e1u(1:jpim1,2:ispongearea) 
    239288    
    240             spe2vr(:      ,2:ispongearea-1) = 0.5 * ( zlocalviscsponge(:,2:ispongearea-1)       & 
    241                &                            +         zlocalviscsponge(:,3:ispongearea  )     ) & 
     289            spe2vr(:      ,2:ispongearea-1) = visc_tra                                     & 
     290               &                            * 0.5 * (  ztabramp(:      ,2:ispongearea-1)   & 
     291               &                                     + ztabramp(:      ,3:ispongearea  ) ) & 
    242292               &                            * e1v(:,2:ispongearea-1) / e2v(:,2:ispongearea-1) 
    243293         ENDIF 
    244294 
    245295         IF( (nbondj == 1) .OR. (nbondj == 2) ) THEN 
    246             DO jj = ilcj+1,nlcj-1 
    247                zlocalviscsponge(:,jj) = visc_tra * (jj - (ilcj+1) ) * z1spongearea 
    248             ENDDO 
    249             spe1ur(1:jpim1,ilcj+1:nlcj-1) = 0.5 * ( zlocalviscsponge(1:jpim1,ilcj+1:nlcj-1)   & 
    250                &                          +         zlocalviscsponge(2:jpi  ,ilcj+1:nlcj-1) ) & 
     296            spe1ur(1:jpim1,ilcj+1:nlcj-1) = visc_tra                                   & 
     297               &                          * 0.5 * (  ztabramp(1:jpim1,ilcj+1:nlcj-1)   & 
     298               &                                   + ztabramp(2:jpi  ,ilcj+1:nlcj-1) ) & 
    251299               &                                * e2u(1:jpim1,ilcj+1:nlcj-1) / e1u(1:jpim1,ilcj+1:nlcj-1) 
    252             spe2vr(:      ,ilcj+1:nlcj-2) = 0.5 * ( zlocalviscsponge(:,ilcj+1:nlcj-2      )   & 
    253                &                          +         zlocalviscsponge(:,ilcj+2:nlcj-1)     )   & 
     300 
     301            spe2vr(:      ,ilcj+1:nlcj-2) = visc_tra                                   & 
     302               &                          * 0.5 * (  ztabramp(:      ,ilcj+1:nlcj-2)   & 
     303               &                                   + ztabramp(:      ,ilcj+2:nlcj-1) ) & 
    254304               &                                * e1v(:,ilcj+1:nlcj-2) / e2v(:,ilcj+1:nlcj-2) 
    255305         ENDIF 
    … …  
    259309      ! Dynamics 
    260310      IF( .NOT. spongedoneU ) THEN 
    261          zlocalviscsponge(:,:) = 0. 
    262311         spe1ur2(:,:) = 0. 
    263312         spe2vr2(:,:) = 0. 
    264313 
    265314         IF( (nbondi == -1) .OR. (nbondi == 2) ) THEN 
    266             DO ji = 2, ispongearea 
    267                zlocalviscsponge(ji,:) = visc_dyn * ( ispongearea-ji ) * z1spongearea 
    268             ENDDO 
    269             spe1ur2(2:ispongearea-1,:      ) = 0.5 * ( zlocalviscsponge(2:ispongearea-1,:      ) & 
    270                                              &     +   zlocalviscsponge(3:ispongearea,:    ) ) 
    271             spe2vr2(2:ispongearea  ,1:jpjm1) = 0.5 * ( zlocalviscsponge(2:ispongearea  ,1:jpjm1) & 
    272                                              &     +   zlocalviscsponge(2:ispongearea,2:jpj) )  
     315            spe1ur2(2:ispongearea-1,:      ) = visc_dyn                                   & 
     316               &                             * 0.5 * (  ztabramp(2:ispongearea-1,:      ) & 
     317               &                                      + ztabramp(3:ispongearea  ,:      ) ) 
     318            spe2vr2(2:ispongearea  ,1:jpjm1) = visc_dyn                                   & 
     319               &                             * 0.5 * (  ztabramp(2:ispongearea  ,1:jpjm1) & 
     320               &                                      + ztabramp(2:ispongearea  ,2:jpj  ) )  
    273321         ENDIF 
    274322 
    275323         IF( (nbondi == 1) .OR. (nbondi == 2) ) THEN 
    276             DO ji = ilci+1,nlci-1 
    277                zlocalviscsponge(ji,:) = visc_dyn * (ji - (ilci+1) ) * z1spongearea 
    278             ENDDO 
    279             spe1ur2(ilci+1:nlci-2,:      ) = 0.5 * (  zlocalviscsponge(ilci+1:nlci-2,:) & 
    280                                            &        + zlocalviscsponge(ilci+2:nlci-1,:) )   
    281             spe2vr2(ilci+1:nlci-1,1:jpjm1) = 0.5 * (  zlocalviscsponge(ilci+1:nlci-1,1:jpjm1) & 
    282                                            &        + zlocalviscsponge(ilci+1:nlci-1,2:jpj  )  )  
     324            spe1ur2(ilci+1:nlci-2  ,:      ) = visc_dyn                                   & 
     325               &                             * 0.5 * (  ztabramp(ilci+1:nlci-2, :       ) & 
     326               &                                      + ztabramp(ilci+2:nlci-1, :       ) )                       
     327            spe2vr2(ilci+1:nlci-1  ,1:jpjm1) = visc_dyn                                   & 
     328               &                             * 0.5 * (  ztabramp(ilci+1:nlci-1,1:jpjm1  ) & 
     329               &                                      + ztabramp(ilci+1:nlci-1,2:jpj    ) )  
    283330         ENDIF 
    284331 
    285332         IF( (nbondj == -1) .OR. (nbondj == 2) ) THEN 
    286             DO jj = 2, ispongearea 
    287                zlocalviscsponge(:,jj) = visc_dyn * ( ispongearea-jj ) * z1spongearea 
    288             ENDDO 
    289             spe1ur2(1:jpim1,2:ispongearea  ) = 0.5 * ( zlocalviscsponge(1:jpim1,2:ispongearea) & 
    290                                              &      + zlocalviscsponge(2:jpi,2:ispongearea) )  
    291             spe2vr2(:      ,2:ispongearea-1) = 0.5 * ( zlocalviscsponge(:,2:ispongearea-1)     & 
    292                                              &      + zlocalviscsponge(:,3:ispongearea)     ) 
     333            spe1ur2(1:jpim1,2:ispongearea  ) = visc_dyn                                   &   
     334               &                             * 0.5 * (  ztabramp(1:jpim1,2:ispongearea  ) & 
     335               &                                      + ztabramp(2:jpi  ,2:ispongearea  ) )  
     336            spe2vr2(:      ,2:ispongearea-1) = visc_dyn                                   & 
     337               &                             * 0.5 * (  ztabramp(:      ,2:ispongearea-1) & 
     338               &                                      + ztabramp(:      ,3:ispongearea  ) ) 
    293339         ENDIF 
    294340 
    295341         IF( (nbondj == 1) .OR. (nbondj == 2) ) THEN 
    296             DO jj = ilcj+1,nlcj-1 
    297                zlocalviscsponge(:,jj) = visc_dyn * (jj - (ilcj+1) ) * z1spongearea 
    298             ENDDO 
    299             spe1ur2(1:jpim1,ilcj+1:nlcj-1) = 0.5 * ( zlocalviscsponge(1:jpim1,ilcj+1:nlcj-1) & 
    300                                            &         + zlocalviscsponge(2:jpi,ilcj+1:nlcj-1) )  
    301             spe2vr2(:      ,ilcj+1:nlcj-2) = 0.5 * ( zlocalviscsponge(:,ilcj+1:nlcj-2      ) & 
    302                                            &         + zlocalviscsponge(:,ilcj+2:nlcj-1)     ) 
     342            spe1ur2(1:jpim1,ilcj+1:nlcj-1  ) = visc_dyn                                   & 
     343               &                             * 0.5 * (  ztabramp(1:jpim1,ilcj+1:nlcj-1  ) & 
     344               &                                      + ztabramp(2:jpi  ,ilcj+1:nlcj-1  ) )  
     345            spe2vr2(:      ,ilcj+1:nlcj-2  ) = visc_dyn                                   & 
     346               &                             * 0.5 * (  ztabramp(:      ,ilcj+1:nlcj-2  ) & 
     347               &                                      + ztabramp(:      ,ilcj+2:nlcj-1  ) ) 
    303348         ENDIF 
    304349         spongedoneU = .TRUE. 
    … …  
    306351      ENDIF 
    307352      ! 
    308       CALL wrk_dealloc( jpi, jpj, zlocalviscsponge ) 
     353      IF (.NOT.ll_spdone) CALL wrk_dealloc( jpi, jpj, ztabramp ) 
    309354      ! 
    310355#endif 

  • branches/2013/dev_UKMO_2013/NEMOGCM/NEMO/OPA_SRC/ASM/asminc.F90

    r3785 r4171  
    682682      ! used to prevent the applied increments taking the temperature below the local freezing point  
    683683 
    684 #if defined key_cice  
    685         fzptnz(:,:,:) = -1.8_wp 
    686 #else  
    687         DO jk = 1, jpk 
    688            DO jj = 1, jpj 
    689               DO ji = 1, jpk 
    690                  fzptnz (ji,jj,jk) = ( -0.0575_wp + 1.710523e-3_wp * SQRT( tsn(ji,jj,jk,jp_sal) )                   &  
    691                                                   - 2.154996e-4_wp *       tsn(ji,jj,jk,jp_sal)   ) * tsn(ji,jj,jk,jp_sal)  &  
    692                                                   - 7.53e-4_wp * fsdepw(ji,jj,jk)       ! (pressure in dbar)  
    693               END DO 
    694            END DO 
    695         END DO 
    696 #endif  
     684      DO jk=1, jpkm1 
     685         fzptnz (:,:,jk) = tfreez( tsn(:,:,jk,jp_sal), fsdept(:,:,jk) ) 
     686      ENDDO 
    697687 
    698688      IF ( ln_asmiau ) THEN 
    … …  
    10261016#endif 
    10271017 
    1028 #if defined key_cice 
     1018#if defined key_cice && defined key_asminc 
    10291019            ! Sea-ice : CICE case. Pass ice increment tendency into CICE 
    10301020            ndaice_da(:,:) = seaice_bkginc(:,:) * zincwgt / rdt 
    … …  
    10371027         ELSE 
    10381028 
    1039 #if defined key_cice 
     1029#if defined key_cice && defined key_asminc 
    10401030            ! Sea-ice : CICE case. Zero ice increment tendency into CICE 
    10411031            ndaice_da(:,:) = 0.0_wp 
    … …  
    10811071#endif 
    10821072  
    1083 #if defined key_cice 
    1084             ! Sea-ice : CICE case. Pass ice increment tendency into CICE - is this correct? 
     1073#if defined key_cice && defined key_asminc 
     1074            ! Sea-ice : CICE case. Pass ice increment tendency into CICE 
    10851075           ndaice_da(:,:) = seaice_bkginc(:,:) / rdt 
    10861076#endif 
    … …  
    10911081         ELSE 
    10921082 
    1093 #if defined key_cice 
     1083#if defined key_cice && defined key_asminc 
    10941084            ! Sea-ice : CICE case. Zero ice increment tendency into CICE  
    10951085            ndaice_da(:,:) = 0.0_wp 

  • branches/2013/dev_UKMO_2013/NEMOGCM/NEMO/OPA_SRC/DIA/diahsb.F90

    r3625 r4171  
    2121   USE bdy_par         ! (for lk_bdy) 
    2222   USE timing          ! preformance summary 
     23   USE lib_fortran 
     24   USE sbcrnf 
    2325 
    2426   IMPLICIT NONE 
    … …  
    3335   REAL(dp)                                ::   surf_tot   , vol_tot             ! 
    3436   REAL(dp)                                ::   frc_t      , frc_s     , frc_v   ! global forcing trends 
     37   REAL(dp)                                ::   frc_wn_t      , frc_wn_s ! global forcing trends 
    3538   REAL(dp)                                ::   fact1                            ! conversion factors 
    3639   REAL(dp)                                ::   fact21    , fact22               !     -         - 
    … …  
    3841   REAL(dp), DIMENSION(:,:)  , ALLOCATABLE ::   surf      , ssh_ini              ! 
    3942   REAL(dp), DIMENSION(:,:,:), ALLOCATABLE ::   hc_loc_ini, sc_loc_ini, e3t_ini  ! 
     43   REAL(dp), DIMENSION(:,:)  , ALLOCATABLE ::   ssh_hc_loc_ini, ssh_sc_loc_ini 
    4044 
    4145   !! * Substitutions 
    … …  
    6771      INTEGER    ::   jk                          ! dummy loop indice 
    6872      REAL(dp)   ::   zdiff_hc    , zdiff_sc      ! heat and salt content variations 
     73      REAL(dp)   ::   zdiff_hc1   , zdiff_sc1     ! heat and salt content variations of ssh 
    6974      REAL(dp)   ::   zdiff_v1    , zdiff_v2      ! volume variation 
     75      REAL(dp)   ::   zerr_hc1    , zerr_sc1      ! Non conservation due to free surface 
    7076      REAL(dp)   ::   z1_rau0                     ! local scalars 
    7177      REAL(dp)   ::   zdeltat                     !    -     - 
    7278      REAL(dp)   ::   z_frc_trd_t , z_frc_trd_s   !    -     - 
    7379      REAL(dp)   ::   z_frc_trd_v                 !    -     - 
     80      REAL(dp)   ::   z_wn_trd_t , z_wn_trd_s   !    -     - 
     81      REAL(dp)   ::   z_ssh_hc , z_ssh_sc   !    -     - 
    7482      !!--------------------------------------------------------------------------- 
    7583      IF( nn_timing == 1 )   CALL timing_start('dia_hsb') 
    … …  
    7987      ! ------------------------- ! 
    8088      z1_rau0 = 1.e0 / rau0 
    81       z_frc_trd_v = z1_rau0 * SUM( - ( emp(:,:) - rnf(:,:) ) * surf(:,:) )     ! volume fluxes 
    82       z_frc_trd_t =           SUM( sbc_tsc(:,:,jp_tem) * surf(:,:) )     ! heat fluxes 
    83       z_frc_trd_s =           SUM( sbc_tsc(:,:,jp_sal) * surf(:,:) )     ! salt fluxes 
     89      z_frc_trd_v = z1_rau0 * glob_sum( - ( emp(:,:) - rnf(:,:) ) * surf(:,:) )     ! volume fluxes 
     90      z_frc_trd_t =           glob_sum( sbc_tsc(:,:,jp_tem) * surf(:,:) )     ! heat fluxes 
     91      z_frc_trd_s =           glob_sum( sbc_tsc(:,:,jp_sal) * surf(:,:) )     ! salt fluxes 
     92      ! Add runoff heat & salt input 
     93      IF( ln_rnf    )   z_frc_trd_t = z_frc_trd_t + glob_sum( rnf_tsc(:,:,jp_tem) * surf(:,:) ) 
     94      IF( ln_rnf_sal)   z_frc_trd_s = z_frc_trd_s + glob_sum( rnf_tsc(:,:,jp_sal) * surf(:,:) ) 
    8495      ! Add penetrative solar radiation 
    85       IF( ln_traqsr )   z_frc_trd_t = z_frc_trd_t + r1_rau0_rcp * SUM( qsr     (:,:) * surf(:,:) ) 
     96      IF( ln_traqsr )   z_frc_trd_t = z_frc_trd_t + r1_rau0_rcp * glob_sum( qsr     (:,:) * surf(:,:) ) 
    8697      ! Add geothermal heat flux 
    87       IF( ln_trabbc )   z_frc_trd_t = z_frc_trd_t + r1_rau0_rcp * SUM( qgh_trd0(:,:) * surf(:,:) ) 
    88       IF( lk_mpp ) THEN 
    89          CALL mpp_sum( z_frc_trd_v ) 
    90          CALL mpp_sum( z_frc_trd_t ) 
    91       ENDIF 
     98      IF( ln_trabbc )   z_frc_trd_t = z_frc_trd_t +  glob_sum( qgh_trd0(:,:) * surf(:,:) ) 
     99      IF( .NOT. lk_vvl ) THEN 
     100         z_wn_trd_t = - glob_sum( surf(:,:) * wn(:,:,1) * tsb(:,:,1,jp_tem) ) 
     101         z_wn_trd_s = - glob_sum( surf(:,:) * wn(:,:,1) * tsb(:,:,1,jp_sal) ) 
     102      ENDIF 
     103 
    92104      frc_v = frc_v + z_frc_trd_v * rdt 
    93105      frc_t = frc_t + z_frc_trd_t * rdt 
    94106      frc_s = frc_s + z_frc_trd_s * rdt 
     107      !                                          ! Advection flux through fixed surface (z=0) 
     108      IF( .NOT. lk_vvl ) THEN 
     109         frc_wn_t = frc_wn_t + z_wn_trd_t * rdt 
     110         frc_wn_s = frc_wn_s + z_wn_trd_s * rdt 
     111      ENDIF 
    95112 
    96113      ! ----------------------- ! 
    … …  
    100117      zdiff_hc = 0.d0 
    101118      zdiff_sc = 0.d0 
     119 
    102120      ! volume variation (calculated with ssh) 
    103       zdiff_v1 = SUM( surf(:,:) * tmask(:,:,1) * ( sshn(:,:) - ssh_ini(:,:) ) ) 
     121      zdiff_v1 = glob_sum( surf(:,:) * ( sshn(:,:) - ssh_ini(:,:) ) ) 
     122 
     123      ! heat & salt content variation (associated with ssh) 
     124      IF( .NOT. lk_vvl ) THEN 
     125         z_ssh_hc = glob_sum( surf(:,:) * ( tsn(:,:,1,jp_tem) * sshn(:,:) - ssh_hc_loc_ini(:,:) ) ) 
     126         z_ssh_sc = glob_sum( surf(:,:) * ( tsn(:,:,1,jp_sal) * sshn(:,:) - ssh_sc_loc_ini(:,:) ) ) 
     127      ENDIF 
     128 
    104129      DO jk = 1, jpkm1 
    105          ! volume variation (calculated with scale factors) 
    106          zdiff_v2 = zdiff_v2 + SUM( surf(:,:) * tmask(:,:,jk)   & 
     130        ! volume variation (calculated with scale factors) 
     131         zdiff_v2 = zdiff_v2 + glob_sum( surf(:,:) * tmask(:,:,jk)   & 
    107132            &                       * ( fse3t_n(:,:,jk)         & 
    108133            &                           - e3t_ini(:,:,jk) ) ) 
    109134         ! heat content variation 
    110          zdiff_hc = zdiff_hc + SUM( surf(:,:) * tmask(:,:,jk)          & 
     135         zdiff_hc = zdiff_hc + glob_sum( surf(:,:) * tmask(:,:,jk)          & 
    111136            &                       * ( fse3t_n(:,:,jk) * tsn(:,:,jk,jp_tem)   & 
    112137            &                           - hc_loc_ini(:,:,jk) ) ) 
    113138         ! salt content variation 
    114          zdiff_sc = zdiff_sc + SUM( surf(:,:) * tmask(:,:,jk)          & 
     139         zdiff_sc = zdiff_sc + glob_sum( surf(:,:) * tmask(:,:,jk)          & 
    115140            &                       * ( fse3t_n(:,:,jk) * tsn(:,:,jk,jp_sal)   & 
    116141            &                           - sc_loc_ini(:,:,jk) ) ) 
    117142      ENDDO 
    118143 
    119       IF( lk_mpp ) THEN 
    120          CALL mpp_sum( zdiff_hc ) 
    121          CALL mpp_sum( zdiff_sc ) 
    122          CALL mpp_sum( zdiff_v1 ) 
    123          CALL mpp_sum( zdiff_v2 ) 
    124       ENDIF 
    125  
    126144      ! Substract forcing from heat content, salt content and volume variations 
    127145      zdiff_v1 = zdiff_v1 - frc_v 
    128       zdiff_v2 = zdiff_v2 - frc_v 
     146      IF( lk_vvl )   zdiff_v2 = zdiff_v2 - frc_v 
    129147      zdiff_hc = zdiff_hc - frc_t 
    130148      zdiff_sc = zdiff_sc - frc_s 
     149      IF( .NOT. lk_vvl ) THEN 
     150         zdiff_hc1 = zdiff_hc + z_ssh_hc  
     151         zdiff_sc1 = zdiff_sc + z_ssh_sc 
     152         zerr_hc1  = z_ssh_hc - frc_wn_t 
     153         zerr_sc1  = z_ssh_sc - frc_wn_s 
     154      ENDIF 
    131155       
    132156      ! ----------------------- ! 
    … …  
    134158      ! ----------------------- ! 
    135159      zdeltat  = 1.e0 / ( ( kt - nit000 + 1 ) * rdt ) 
    136       WRITE(numhsb , 9020) kt , zdiff_hc / vol_tot , zdiff_hc * fact1  * zdeltat,                                & 
    137          &                      zdiff_sc / vol_tot , zdiff_sc * fact21 * zdeltat, zdiff_sc * fact22 * zdeltat,   & 
    138          &                      zdiff_v1           , zdiff_v1 * fact31 * zdeltat, zdiff_v1 * fact32 * zdeltat,   & 
    139          &                      zdiff_v2           , zdiff_v2 * fact31 * zdeltat, zdiff_v2 * fact32 * zdeltat 
     160      IF( lk_vvl ) THEN 
     161         WRITE(numhsb , 9020) kt , zdiff_hc / vol_tot , zdiff_hc * fact1  * zdeltat,                                & 
     162            &                      zdiff_sc / vol_tot , zdiff_sc * fact21 * zdeltat, zdiff_sc * fact22 * zdeltat,   & 
     163            &                      zdiff_v1           , zdiff_v1 * fact31 * zdeltat, zdiff_v1 * fact32 * zdeltat,   & 
     164            &                      zdiff_v2           , zdiff_v2 * fact31 * zdeltat, zdiff_v2 * fact32 * zdeltat 
     165      ELSE 
     166         WRITE(numhsb , 9030) kt , zdiff_hc1 / vol_tot , zdiff_hc1 * fact1  * zdeltat,                                & 
     167            &                      zdiff_sc1 / vol_tot , zdiff_sc1 * fact21 * zdeltat, zdiff_sc1 * fact22 * zdeltat,   & 
     168            &                      zdiff_v1            , zdiff_v1  * fact31 * zdeltat, zdiff_v1  * fact32 * zdeltat,   & 
     169            &                      zerr_hc1 / vol_tot  , zerr_sc1 / vol_tot 
     170      ENDIF 
    140171 
    141172      IF ( kt == nitend ) CLOSE( numhsb ) 
    … …  
    144175 
    1451769020  FORMAT(I5,11D15.7) 
     1779030  FORMAT(I5,10D15.7) 
    146178      ! 
    147179   END SUBROUTINE dia_hsb 
    … …  
    179211 
    180212      IF( .NOT. ln_diahsb )   RETURN 
     213      IF( .NOT. lk_mpp_rep ) & 
     214        CALL ctl_stop (' Your global mpp_sum if performed in single precision - 64 bits -', & 
     215             &         ' whereas the global sum to be precise must be done in double precision ',& 
     216             &         ' please add key_mpp_rep') 
    181217 
    182218      ! ------------------- ! 
    183219      ! 1 - Allocate memory ! 
    184220      ! ------------------- ! 
    185       ALLOCATE( hc_loc_ini(jpi,jpj,jpk), STAT=ierror ) 
     221      ALLOCATE( hc_loc_ini(jpi,jpj,jpk), sc_loc_ini(jpi,jpj,jpk), & 
     222         &      ssh_hc_loc_ini(jpi,jpj), ssh_sc_loc_ini(jpi,jpj), & 
     223         &      e3t_ini(jpi,jpj,jpk)                            , & 
     224         &      surf(jpi,jpj),  ssh_ini(jpi,jpj), STAT=ierror ) 
    186225      IF( ierror > 0 ) THEN 
    187226         CALL ctl_stop( 'dia_hsb: unable to allocate hc_loc_ini' )   ;   RETURN 
    188       ENDIF 
    189       ALLOCATE( sc_loc_ini(jpi,jpj,jpk), STAT=ierror ) 
    190       IF( ierror > 0 ) THEN 
    191          CALL ctl_stop( 'dia_hsb: unable to allocate sc_loc_ini' )   ;   RETURN 
    192       ENDIF 
    193       ALLOCATE( e3t_ini(jpi,jpj,jpk)   , STAT=ierror ) 
    194       IF( ierror > 0 ) THEN 
    195          CALL ctl_stop( 'dia_hsb: unable to allocate e3t_ini' )      ;   RETURN 
    196       ENDIF 
    197       ALLOCATE( surf(jpi,jpj)          , STAT=ierror ) 
    198       IF( ierror > 0 ) THEN 
    199          CALL ctl_stop( 'dia_hsb: unable to allocate surf' )         ;   RETURN 
    200       ENDIF 
    201       ALLOCATE( ssh_ini(jpi,jpj)       , STAT=ierror ) 
    202       IF( ierror > 0 ) THEN 
    203          CALL ctl_stop( 'dia_hsb: unable to allocate ssh_ini' )      ;   RETURN 
    204227      ENDIF 
    205228 
    … …  
    214237      cl_name    = 'heat_salt_volume_budgets.txt'                         ! name of output file 
    215238      surf(:,:) = e1t(:,:) * e2t(:,:) * tmask(:,:,1) * tmask_i(:,:)      ! masked surface grid cell area 
    216       surf_tot  = SUM( surf(:,:) )                                       ! total ocean surface area 
     239      surf_tot  = glob_sum( surf(:,:) )                                       ! total ocean surface area 
    217240      vol_tot   = 0.d0                                                   ! total ocean volume 
    218241      DO jk = 1, jpkm1 
    219          vol_tot  = vol_tot + SUM( surf(:,:) * tmask(:,:,jk)     & 
    220             &                      * fse3t_n(:,:,jk)         ) 
     242         vol_tot  = vol_tot + glob_sum( surf(:,:) * tmask(:,:,jk)     & 
     243            &                         * fse3t_n(:,:,jk)         ) 
    221244      END DO 
    222       IF( lk_mpp ) THEN  
    223          CALL mpp_sum( vol_tot ) 
    224          CALL mpp_sum( surf_tot ) 
    225       ENDIF 
    226245 
    227246      CALL ctl_opn( numhsb , cl_name , 'UNKNOWN' , 'FORMATTED' , 'SEQUENTIAL' , 1 , numout , lwp , 1 ) 
    228       !                   12345678901234567890123456789012345678901234567890123456789012345678901234567890 -> 80 
    229       WRITE( numhsb, 9010 ) "kt   |     heat content budget     |            salt content budget             ",   & 
    230          !                                                   123456789012345678901234567890123456789012345 -> 45 
    231          &                                                  "|            volume budget (ssh)             ",   & 
    232          !                                                   678901234567890123456789012345678901234567890 -> 45 
    233          &                                                  "|            volume budget (e3t)             " 
    234       WRITE( numhsb, 9010 ) "     |      [C]         [W/m2]     |     [psu]        [mmm/s]          [SV]     ",   & 
    235          &                                                  "|     [m3]         [mmm/s]          [SV]     ",   & 
    236          &                                                  "|     [m3]         [mmm/s]          [SV]     " 
    237  
     247      IF( lk_vvl ) THEN 
     248         !                   12345678901234567890123456789012345678901234567890123456789012345678901234567890 -> 80 
     249         WRITE( numhsb, 9010 ) "kt   |     heat content budget     |            salt content budget             ",   & 
     250            !                                                   123456789012345678901234567890123456789012345 -> 45 
     251            &                                                  "|            volume budget (ssh)             ",   & 
     252            !                                                   678901234567890123456789012345678901234567890 -> 45 
     253            &                                                  "|            volume budget (e3t)             " 
     254         WRITE( numhsb, 9010 ) "     |      [C]         [W/m2]     |     [psu]        [mmm/s]          [SV]     ",   & 
     255            &                                                  "|     [m3]         [mmm/s]          [SV]     ",   & 
     256            &                                                  "|     [m3]         [mmm/s]          [SV]     " 
     257      ELSE 
     258         !                   12345678901234567890123456789012345678901234567890123456789012345678901234567890 -> 80 
     259         WRITE( numhsb, 9011 ) "kt   |     heat content budget     |            salt content budget             ",   & 
     260            !                                                   123456789012345678901234567890123456789012345 -> 45 
     261            &                                                  "|            volume budget (ssh)             ",   & 
     262            !                                                   678901234567890123456789012345678901234567890 -> 45 
     263            &                                                  "|  Non conservation due to free surface      " 
     264         WRITE( numhsb, 9011 ) "     |      [C]         [W/m2]     |     [psu]        [mmm/s]          [SV]     ",   & 
     265            &                                                  "|     [m3]         [mmm/s]          [SV]     ",   & 
     266            &                                                  "|  [heat - C]     [salt - psu]                " 
     267      ENDIF 
    238268      ! --------------- ! 
    239269      ! 3 - Conversions ! (factors will be multiplied by duration afterwards) 
    … …  
    261291      frc_t = 0.d0                                           ! heat content   -    -   -    -    
    262292      frc_s = 0.d0                                           ! salt content   -    -   -    -          
     293      IF( .NOT. lk_vvl ) THEN 
     294         ssh_hc_loc_ini(:,:) = tsn(:,:,1,jp_tem) * ssh_ini(:,:)   ! initial heat content associated with ssh 
     295         ssh_sc_loc_ini(:,:) = tsn(:,:,1,jp_sal) * ssh_ini(:,:)   ! initial salt content associated with ssh 
     296         frc_wn_t = 0.d0 
     297         frc_wn_s = 0.d0 
     298      ENDIF 
    263299      ! 
    2643009010  FORMAT(A80,A45,A45) 
     3019011  FORMAT(A80,A45,A45) 
    265302      ! 
    266303   END SUBROUTINE dia_hsb_init 

  • branches/2013/dev_UKMO_2013/NEMOGCM/NEMO/OPA_SRC/DOM/closea.F90

    r3632 r4171  
    108108            ncsi1(2)   =  97  ;  ncsj1(2)   = 107 
    109109            ncsi2(2)   = 103  ;  ncsj2(2)   = 111 
    110             ncsir(2,1) = 110  ;  ncsjr(2,1) = 111 
    111             !                                            ! Black Sea 1 : west part of the Black Sea  
    112             ncsnr(3)   = 1    ; ncstt(3)   =   2            !            (ie west of the cyclic b.c.) 
    113             ncsi1(3)   = 174  ; ncsj1(3)   = 107            ! put in Med Sea 
    114             ncsi2(3)   = 181  ; ncsj2(3)   = 112 
    115             ncsir(3,1) = 171  ; ncsjr(3,1) = 106  
    116             !                                            ! Black Sea 2 : est part of the Black Sea  
    117             ncsnr(4)   =   1  ;  ncstt(4)   =   2           !               (ie est of the cyclic b.c.) 
    118             ncsi1(4)   =   2  ;  ncsj1(4)   = 107           ! put in Med Sea 
    119             ncsi2(4)   =   6  ;  ncsj2(4)   = 112 
    120             ncsir(4,1) = 171  ;  ncsjr(4,1) = 106  
     110            ncsir(2,1) = 110  ;  ncsjr(2,1) = 111            
     111            !                                            ! Black Sea (crossed by the cyclic boundary condition) 
     112            ncsnr(3:4) =   4  ;  ncstt(3:4) =   2           ! put in Med Sea (north of Aegean Sea) 
     113            ncsir(3:4,1) = 171;  ncsjr(3:4,1) = 106         ! 
     114            ncsir(3:4,2) = 170;  ncsjr(3:4,2) = 106  
     115            ncsir(3:4,3) = 171;  ncsjr(3:4,3) = 105  
     116            ncsir(3:4,4) = 170;  ncsjr(3:4,4) = 105  
     117            ncsi1(3)   = 174  ;  ncsj1(3)   = 107           ! 1 : west part of the Black Sea       
     118            ncsi2(3)   = 181  ;  ncsj2(3)   = 112           !            (ie west of the cyclic b.c.) 
     119            ncsi1(4)   =   2  ;  ncsj1(4)   = 107           ! 2 : east part of the Black Sea  
     120            ncsi2(4)   =   6  ;  ncsj2(4)   = 112           !           (ie east of the cyclic b.c.) 
     121              
     122           
     123 
    121124            !                                        ! ======================= 
    122125         CASE ( 4 )                                  !  ORCA_R4 configuration 
    … …  
    372375      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   p_rnfmsk   ! river runoff mask (rnfmsk array) 
    373376      ! 
    374       INTEGER  ::   jc, jn      ! dummy loop indices 
    375       INTEGER  ::   ii, ij      ! temporary integer 
     377      INTEGER  ::   jc, jn, ji, jj      ! dummy loop indices 
    376378      !!---------------------------------------------------------------------- 
    377379      ! 
    … …  
    379381         IF( ncstt(jc) >= 1 ) THEN            ! runoff mask set to 1 at closed sea outflows 
    380382             DO jn = 1, 4 
    381                ii = mi0( ncsir(jc,jn) ) 
    382                ij = mj0( ncsjr(jc,jn) ) 
    383                p_rnfmsk(ii,ij) = MAX( p_rnfmsk(ii,ij), 1.0_wp ) 
     383                DO jj =    mj0( ncsjr(jc,jn) ), mj1( ncsjr(jc,jn) ) 
     384                   DO ji = mi0( ncsir(jc,jn) ), mi1( ncsir(jc,jn) ) 
     385                      p_rnfmsk(ji,jj) = MAX( p_rnfmsk(ji,jj), 1.0_wp ) 
     386                   END DO 
     387                END DO 
    384388            END DO  
    385389         ENDIF  

  • branches/2013/dev_UKMO_2013/NEMOGCM/NEMO/OPA_SRC/DOM/daymod.F90

    r3851 r4171  
    238238               nday_year = 1 
    239239               nsec_year = ndt05 
     240               IF( nsec1jan000 >= 2 * (2**30 - nsecd * nyear_len(1) / 2 ) ) THEN   ! test integer 4 max value 
     241                  CALL ctl_stop( 'The number of seconds between Jan. 1st 00h of nit000 year and Jan. 1st 00h ',   & 
     242                     &           'of the current year is exceeding the INTEGER 4 max VALUE: 2^31-1 -> 68.09 years in seconds', & 
     243                     & 'You must do a restart at higher frequency (or remove this STOP and recompile everything in I8)' ) 
     244               ENDIF 
    240245               nsec1jan000 = nsec1jan000 + nsecd * nyear_len(1) 
    241246               IF( nleapy == 1 )   CALL day_mth 

  • branches/2013/dev_UKMO_2013/NEMOGCM/NEMO/OPA_SRC/DOM/domzgr.F90

    r3971 r4171  
    11021102      INTEGER  ::   iip1, ijp1, iim1, ijm1   ! temporary integers 
    11031103      REAL(wp) ::   zrmax, ztaper   ! temporary scalars 
    1104       REAL(wp) ::   zrfact   ! temporary scalars 
    1105       REAL(wp), POINTER, DIMENSION(:,:  ) :: ztmpi1, ztmpi2, ztmpj1, ztmpj2 
    1106  
    1107       ! 
    1108       REAL(wp), POINTER, DIMENSION(:,:  ) :: zenv, zri, zrj, zhbat 
     1104      ! 
     1105      REAL(wp), POINTER, DIMENSION(:,:  ) :: zenv, ztmp, zmsk, zri, zrj, zhbat 
    11091106 
    11101107      NAMELIST/namzgr_sco/ln_s_sh94, ln_s_sf12, ln_sigcrit, rn_sbot_min, rn_sbot_max, rn_hc, rn_rmax,rn_theta, & 
    … …  
    11141111      IF( nn_timing == 1 )  CALL timing_start('zgr_sco') 
    11151112      ! 
    1116       CALL wrk_alloc( jpi, jpj,      ztmpi1, ztmpi2, ztmpj1, ztmpj2         ) 
    1117       CALL wrk_alloc( jpi, jpj,      zenv, zri, zrj, zhbat     ) 
    1118      ! 
     1113      CALL wrk_alloc( jpi, jpj,      zenv, ztmp, zmsk, zri, zrj, zhbat                           ) 
     1114      ! 
    11191115      REWIND( numnam )                       ! Read Namelist namzgr_sco : sigma-stretching parameters 
    11201116      READ  ( numnam, namzgr_sco ) 
    … …  
    11631159      !                                        ! ============================= 
    11641160      ! use r-value to create hybrid coordinates 
    1165 !     DO jj = 1, jpj 
    1166 !        DO ji = 1, jpi 
    1167 !           zenv(ji,jj) = MAX( bathy(ji,jj), 0._wp ) 
    1168 !        END DO 
    1169 !     END DO 
    1170 !     CALL lbc_lnk( zenv, 'T', 1._wp ) 
    1171       zenv(:,:) = bathy(:,:) 
     1161      DO jj = 1, jpj 
     1162         DO ji = 1, jpi 
     1163            zenv(ji,jj) = MAX( bathy(ji,jj), rn_sbot_min ) 
     1164         END DO 
     1165      END DO 
    11721166      !  
    11731167      ! Smooth the bathymetry (if required) 
    … …  
    11771171      jl = 0 
    11781172      zrmax = 1._wp 
    1179       !      
    1180       ! set scaling factor used in reducing vertical gradients 
    1181       zrfact = ( 1._wp - rn_rmax ) / ( 1._wp + rn_rmax )  
    1182       ! 
    1183       ! initialise temporary evelope depth arrays 
    1184       ztmpi1(:,:) = zenv(:,:) 
    1185       ztmpi2(:,:) = zenv(:,:) 
    1186       ztmpj1(:,:) = zenv(:,:) 
    1187       ztmpj2(:,:) = zenv(:,:) 
    1188       ! 
    1189       ! initialise temporary r-value arrays 
    1190       zri(:,:) = 1._wp 
    1191       zrj(:,:) = 1._wp 
    1192       !                                                            ! ================ ! 
    1193       DO WHILE( jl <= 10000 .AND. ( zrmax - rn_rmax ) > 1.e-8_wp ) !  Iterative loop  ! 
    1194          !                                                         ! ================ ! 
     1173      !                                                     ! ================ ! 
     1174      DO WHILE( jl <= 10000 .AND. zrmax > rn_rmax )         !  Iterative loop  ! 
     1175         !                                                  ! ================ ! 
    11951176         jl = jl + 1 
    11961177         zrmax = 0._wp 
    1197          ! we set zrmax from previous r-values (zri abd zrj) first 
    1198          ! if set after current r-value calculation (as previously) 
    1199          ! we could exit DO WHILE prematurely before checking r-value 
    1200          ! of current zenv 
    1201          DO jj = 1, nlcj 
    1202             DO ji = 1, nlci 
    1203                zrmax = MAX( zrmax, ABS(zri(ji,jj)), ABS(zrj(ji,jj)) ) 
    1204             END DO 
    1205          END DO 
    1206          zri(:,:) = 0._wp 
    1207          zrj(:,:) = 0._wp 
     1178         zmsk(:,:) = 0._wp 
    12081179         DO jj = 1, nlcj 
    12091180            DO ji = 1, nlci 
    12101181               iip1 = MIN( ji+1, nlci )      ! force zri = 0 on last line (ji=ncli+1 to jpi) 
    12111182               ijp1 = MIN( jj+1, nlcj )      ! force zrj = 0 on last raw  (jj=nclj+1 to jpj) 
    1212                IF( (zenv(ji,jj) > 0._wp) .AND. (zenv(iip1,jj) > 0._wp)) THEN 
    1213                   zri(ji,jj) = ( zenv(iip1,jj  ) - zenv(ji,jj) ) / ( zenv(iip1,jj  ) + zenv(ji,jj) ) 
    1214                END IF 
    1215                IF( (zenv(ji,jj) > 0._wp) .AND. (zenv(ji,ijp1) > 0._wp)) THEN 
    1216                   zrj(ji,jj) = ( zenv(ji  ,ijp1) - zenv(ji,jj) ) / ( zenv(ji  ,ijp1) + zenv(ji,jj) ) 
    1217                END IF 
    1218                IF( zri(ji,jj) >  rn_rmax )   ztmpi1(ji  ,jj  ) = zenv(iip1,jj  ) * zrfact 
    1219                IF( zri(ji,jj) < -rn_rmax )   ztmpi2(iip1,jj  ) = zenv(ji  ,jj  ) * zrfact  
    1220                IF( zrj(ji,jj) >  rn_rmax )   ztmpj1(ji  ,jj  ) = zenv(ji  ,ijp1) * zrfact 
    1221                IF( zrj(ji,jj) < -rn_rmax )   ztmpj2(ji  ,ijp1) = zenv(ji  ,jj  ) * zrfact 
     1183               zri(ji,jj) = ABS( zenv(iip1,jj  ) - zenv(ji,jj) ) / ( zenv(iip1,jj  ) + zenv(ji,jj) ) 
     1184               zrj(ji,jj) = ABS( zenv(ji  ,ijp1) - zenv(ji,jj) ) / ( zenv(ji  ,ijp1) + zenv(ji,jj) ) 
     1185               zrmax = MAX( zrmax, zri(ji,jj), zrj(ji,jj) ) 
     1186               IF( zri(ji,jj) > rn_rmax )   zmsk(ji  ,jj  ) = 1._wp 
     1187               IF( zri(ji,jj) > rn_rmax )   zmsk(iip1,jj  ) = 1._wp 
     1188               IF( zrj(ji,jj) > rn_rmax )   zmsk(ji  ,jj  ) = 1._wp 
     1189               IF( zrj(ji,jj) > rn_rmax )   zmsk(ji  ,ijp1) = 1._wp 
    12221190            END DO 
    12231191         END DO 
    12241192         IF( lk_mpp )   CALL mpp_max( zrmax )   ! max over the global domain 
     1193         ! lateral boundary condition on zmsk: keep 1 along closed boundary (use of MAX) 
     1194         ztmp(:,:) = zmsk(:,:)   ;   CALL lbc_lnk( zmsk, 'T', 1._wp ) 
     1195         DO jj = 1, nlcj 
     1196            DO ji = 1, nlci 
     1197                zmsk(ji,jj) = MAX( zmsk(ji,jj), ztmp(ji,jj) ) 
     1198            END DO 
     1199         END DO 
    12251200         ! 
    1226          IF(lwp)WRITE(numout,*) 'zgr_sco :   iter= ',jl, ' rmax= ', zrmax 
     1201         IF(lwp)WRITE(numout,*) 'zgr_sco :   iter= ',jl, ' rmax= ', zrmax, ' nb of pt= ', INT( SUM(zmsk(:,:) ) ) 
    12271202         ! 
    12281203         DO jj = 1, nlcj 
    12291204            DO ji = 1, nlci 
    1230                zenv(ji,jj) = MAX(zenv(ji,jj), ztmpi1(ji,jj), ztmpi2(ji,jj), ztmpj1(ji,jj), ztmpj2(ji,jj) ) 
     1205               iip1 = MIN( ji+1, nlci )     ! last  line (ji=nlci) 
     1206               ijp1 = MIN( jj+1, nlcj )     ! last  raw  (jj=nlcj) 
     1207               iim1 = MAX( ji-1,  1  )      ! first line (ji=nlci) 
     1208               ijm1 = MAX( jj-1,  1  )      ! first raw  (jj=nlcj) 
     1209               IF( zmsk(ji,jj) == 1._wp ) THEN 
     1210                  ztmp(ji,jj) =   (                                                                                   & 
     1211             &      zenv(iim1,ijp1)*zmsk(iim1,ijp1) + zenv(ji,ijp1)*zmsk(ji,ijp1) + zenv(iip1,ijp1)*zmsk(iip1,ijp1)   & 
     1212             &    + zenv(iim1,jj  )*zmsk(iim1,jj  ) + zenv(ji,jj  )*    2._wp     + zenv(iip1,jj  )*zmsk(iip1,jj  )   & 
     1213             &    + zenv(iim1,ijm1)*zmsk(iim1,ijm1) + zenv(ji,ijm1)*zmsk(ji,ijm1) + zenv(iip1,ijm1)*zmsk(iip1,ijm1)   & 
     1214             &                    ) / (                                                                               & 
     1215             &                      zmsk(iim1,ijp1) +               zmsk(ji,ijp1) +                 zmsk(iip1,ijp1)   & 
     1216             &    +                 zmsk(iim1,jj  ) +                   2._wp     +                 zmsk(iip1,jj  )   & 
     1217             &    +                 zmsk(iim1,ijm1) +               zmsk(ji,ijm1) +                 zmsk(iip1,ijm1)   & 
     1218             &                        ) 
     1219               ENDIF 
    12311220            END DO 
    12321221         END DO 
    12331222         ! 
    1234          CALL lbc_lnk( zenv, 'T', 1._wp ) 
     1223         DO jj = 1, nlcj 
     1224            DO ji = 1, nlci 
     1225               IF( zmsk(ji,jj) == 1._wp )   zenv(ji,jj) = MAX( ztmp(ji,jj), bathy(ji,jj) ) 
     1226            END DO 
     1227         END DO 
     1228         ! 
     1229         ! Apply lateral boundary condition   CAUTION: keep the value when the lbc field is zero 
     1230         ztmp(:,:) = zenv(:,:)   ;   CALL lbc_lnk( zenv, 'T', 1._wp ) 
     1231         DO jj = 1, nlcj 
     1232            DO ji = 1, nlci 
     1233               IF( zenv(ji,jj) == 0._wp )   zenv(ji,jj) = ztmp(ji,jj) 
     1234            END DO 
     1235         END DO 
    12351236         !                                                  ! ================ ! 
    12361237      END DO                                                !     End loop     ! 
    12371238      !                                                     ! ================ ! 
    12381239      ! 
    1239 !     DO jj = 1, jpj 
    1240 !        DO ji = 1, jpi 
    1241 !           zenv(ji,jj) = MAX( zenv(ji,jj), rn_sbot_min ) ! set all points to avoid undefined scale values 
    1242 !        END DO 
    1243 !     END DO 
     1240      ! Fill ghost rows with appropriate values to avoid undefined e3 values with some mpp decompositions 
     1241      DO ji = nlci+1, jpi  
     1242         zenv(ji,1:nlcj) = zenv(nlci,1:nlcj) 
     1243      END DO 
     1244      ! 
     1245      DO jj = nlcj+1, jpj 
     1246         zenv(:,jj) = zenv(:,nlcj) 
     1247      END DO 
    12441248      ! 
    12451249      ! Envelope bathymetry saved in hbatt 
    12461250      hbatt(:,:) = zenv(:,:)  
    1247  
    12481251      IF( MINVAL( gphit(:,:) ) * MAXVAL( gphit(:,:) ) <= 0._wp ) THEN 
    12491252         CALL ctl_warn( ' s-coordinates are tapered in vicinity of the Equator' ) 
    12501253         DO jj = 1, jpj 
    12511254            DO ji = 1, jpi 
    1252                ztaper = EXP( -(gphit(ji,jj)/8._wp)**2 ) 
     1255               ztaper = EXP( -(gphit(ji,jj)/8._wp)**2._wp ) 
    12531256               hbatt(ji,jj) = rn_sbot_max * ztaper + hbatt(ji,jj) * ( 1._wp - ztaper ) 
    12541257            END DO 
    … …  
    13651368      fsde3w(:,:,:) = gdep3w(:,:,:) 
    13661369      ! 
    1367       where (e3t   (:,:,:).eq.0.0)  e3t(:,:,:) = 1.0 
    1368       where (e3u   (:,:,:).eq.0.0)  e3u(:,:,:) = 1.0 
    1369       where (e3v   (:,:,:).eq.0.0)  e3v(:,:,:) = 1.0 
    1370       where (e3f   (:,:,:).eq.0.0)  e3f(:,:,:) = 1.0 
    1371       where (e3w   (:,:,:).eq.0.0)  e3w(:,:,:) = 1.0 
    1372       where (e3uw  (:,:,:).eq.0.0)  e3uw(:,:,:) = 1.0 
    1373       where (e3vw  (:,:,:).eq.0.0)  e3vw(:,:,:) = 1.0 
     1370      where (e3t   (:,:,:).eq.0.0)  e3t(:,:,:) = 1._wp 
     1371      where (e3u   (:,:,:).eq.0.0)  e3u(:,:,:) = 1._wp 
     1372      where (e3v   (:,:,:).eq.0.0)  e3v(:,:,:) = 1._wp 
     1373      where (e3f   (:,:,:).eq.0.0)  e3f(:,:,:) = 1._wp 
     1374      where (e3w   (:,:,:).eq.0.0)  e3w(:,:,:) = 1._wp 
     1375      where (e3uw  (:,:,:).eq.0.0)  e3uw(:,:,:) = 1._wp 
     1376      where (e3vw  (:,:,:).eq.0.0)  e3vw(:,:,:) = 1._wp 
    13741377 
    13751378#if defined key_agrif 
    … …  
    15191522      END DO 
    15201523      ! 
    1521       CALL wrk_dealloc( jpi, jpj,      zenv, ztmpi1, ztmpi2, ztmpj1, ztmpj2, zri, zrj, zhbat                           )      ! 
     1524      CALL wrk_dealloc( jpi, jpj,      zenv, ztmp, zmsk, zri, zrj, zhbat                           ) 
     1525      ! 
    15221526      IF( nn_timing == 1 )  CALL timing_stop('zgr_sco') 
    15231527      ! 
    … …  
    17481752      ENDDO 
    17491753      ! 
    1750       CALL lbc_lnk(e3t ,'T',1.) ; CALL lbc_lnk(e3u ,'T',1.) 
    1751       CALL lbc_lnk(e3v ,'T',1.) ; CALL lbc_lnk(e3f ,'T',1.) 
    1752       CALL lbc_lnk(e3w ,'T',1.) 
    1753       CALL lbc_lnk(e3uw,'T',1.) ; CALL lbc_lnk(e3vw,'T',1.) 
    1754       ! 
    17551754      !                                               ! ============= 
    17561755 
    … …  
    18491848      !!---------------------------------------------------------------------- 
    18501849      ! 
    1851       pf =   (   TANH( rn_theta * ( -(pk-0.5_wp) / REAL(jpkm1) + rn_thetb )  )   & 
     1850      pf =   (   TANH( rn_theta * ( -(pk-0.5_wp) / REAL(jpkm1,wp) + rn_thetb )  )   & 
    18521851         &     - TANH( rn_thetb * rn_theta                                )  )   & 
    18531852         & * (   COSH( rn_theta                           )                      & 
    … …  
    18751874      ! 
    18761875      IF ( rn_theta == 0 ) then      ! uniform sigma 
    1877          pf1 = - ( pk1 - 0.5_wp ) / REAL( jpkm1 ) 
     1876         pf1 = - ( pk1 - 0.5_wp ) / REAL( jpkm1,wp ) 
    18781877      ELSE                        ! stretched sigma 
    1879          pf1 =   ( 1._wp - pbb ) * ( SINH( rn_theta*(-(pk1-0.5_wp)/REAL(jpkm1)) ) ) / SINH( rn_theta )              & 
    1880             &  + pbb * (  (TANH( rn_theta*( (-(pk1-0.5_wp)/REAL(jpkm1)) + 0.5_wp) ) - TANH( 0.5_wp * rn_theta )  )  & 
     1878         pf1 =   ( 1._wp - pbb ) * ( SINH( rn_theta*(-(pk1-0.5_wp)/REAL(jpkm1,wp)) ) ) / SINH( rn_theta )              & 
     1879            &  + pbb * (  (TANH( rn_theta*( (-(pk1-0.5_wp)/REAL(jpkm1,wp)) + 0.5_wp) ) - TANH( 0.5_wp * rn_theta )  )  & 
    18811880            &        / ( 2._wp * TANH( 0.5_wp * rn_theta ) )  ) 
    18821881      ENDIF 

  • branches/2013/dev_UKMO_2013/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_flt.F90

    r3765 r4171  
    109109      INTEGER  ::   ji, jj, jk   ! dummy loop indices 
    110110      REAL(wp) ::   z2dt, z2dtg, zgcb, zbtd, ztdgu, ztdgv   ! local scalars 
    111       REAL(wp), POINTER, DIMENSION(:,:,:) ::  zub, zvb 
    112111      !!---------------------------------------------------------------------- 
    113112      ! 
    114113      IF( nn_timing == 1 )  CALL timing_start('dyn_spg_flt') 
    115114      ! 
    116       CALL wrk_alloc( jpi,jpj,jpk, zub, zvb ) 
    117115      ! 
    118116      IF( kt == nit000 ) THEN 
    … …  
    213211         DO jk = 1, jpkm1 
    214212            DO ji = 1, jpij 
    215                spgu(ji,1) = spgu(ji,1) + fse3u(ji,1,jk) * ua(ji,1,jk) 
    216                spgv(ji,1) = spgv(ji,1) + fse3v(ji,1,jk) * va(ji,1,jk) 
     213               spgu(ji,1) = spgu(ji,1) + fse3u_a(ji,1,jk) * ua(ji,1,jk) 
     214               spgv(ji,1) = spgv(ji,1) + fse3v_a(ji,1,jk) * va(ji,1,jk) 
    217215            END DO 
    218216         END DO 
    … …  
    221219            DO jj = 2, jpjm1 
    222220               DO ji = 2, jpim1 
    223                   spgu(ji,jj) = spgu(ji,jj) + fse3u(ji,jj,jk) * ua(ji,jj,jk) 
    224                   spgv(ji,jj) = spgv(ji,jj) + fse3v(ji,jj,jk) * va(ji,jj,jk) 
     221                  spgu(ji,jj) = spgu(ji,jj) + fse3u_a(ji,jj,jk) * ua(ji,jj,jk) 
     222                  spgv(ji,jj) = spgv(ji,jj) + fse3v_a(ji,jj,jk) * va(ji,jj,jk) 
    225223               END DO 
    226224            END DO 
    … …  
    360358      IF( lrst_oce ) CALL flt_rst( kt, 'WRITE' ) 
    361359      ! 
    362       CALL wrk_dealloc( jpi,jpj,jpk, zub, zvb ) 
    363360      ! 
    364361      IF( nn_timing == 1 )  CALL timing_stop('dyn_spg_flt') 

  • branches/2013/dev_UKMO_2013/NEMOGCM/NEMO/OPA_SRC/ICB/icb_oce.F90

    r3983 r4171  
    165165      ! 
    166166      icb_alloc = 0 
    167       ALLOCATE( berg_grid                      ,                                               & 
    168          &      berg_grid%calving    (jpi,jpj) , berg_grid%calving_hflx (jpi,jpj)          ,   & 
     167      ALLOCATE( berg_grid%calving    (jpi,jpj) , berg_grid%calving_hflx (jpi,jpj)          ,   & 
    169168         &      berg_grid%stored_heat(jpi,jpj) , berg_grid%floating_melt(jpi,jpj)          ,   & 
    170169         &      berg_grid%maxclass   (jpi,jpj) , berg_grid%stored_ice   (jpi,jpj,nclasses) ,   & 

  • branches/2013/dev_UKMO_2013/NEMOGCM/NEMO/OPA_SRC/LBC/lib_mpp.F90

    r3918 r4171  
    21792179!!gm Remark : this is very time consumming!!! 
    21802180      !                                         ! ------------------------ ! 
    2181             IF( ijpt0 > ijpt1 .OR. iipt0 > iipt1 ) THEN 
     2181        IF(((nbondi .ne. 0) .AND. (ktype .eq. 2)) .OR. ((nbondj .ne. 0) .AND. (ktype .eq. 1))) THEN 
    21822182            ! there is nothing to be migrated 
    2183                lmigr = .FALSE. 
     2183              lmigr = .TRUE. 
    21842184            ELSE 
    2185               lmigr = .TRUE. 
     2185              lmigr = .FALSE. 
    21862186            ENDIF 
    21872187 

  • branches/2013/dev_UKMO_2013/NEMOGCM/NEMO/OPA_SRC/OBS/diaobs.F90

    r3651 r4171  
    2121   USE par_oce 
    2222   USE dom_oce                  ! Ocean space and time domain variables 
     23   USE obs_fbm, ONLY: ln_cl4    ! Class 4 diagnostic switch 
    2324   USE obs_read_prof            ! Reading and allocation of observations (Coriolis) 
    2425   USE obs_read_sla             ! Reading and allocation of SLA observations   
    … …  
    4849   PUBLIC dia_obs_init, &  ! Initialize and read observations 
    4950      &   dia_obs,      &  ! Compute model equivalent to observations 
    50       &   dia_obs_wri      ! Write model equivalent to observations 
     51      &   dia_obs_wri,  &  ! Write model equivalent to observations 
     52      &   dia_obs_dealloc  ! Deallocate dia_obs data 
    5153 
    5254   !! * Shared Module variables 
    … …  
    8082   LOGICAL, PUBLIC :: ln_ssh         !: Logical switch for sea surface height 
    8183   LOGICAL, PUBLIC :: ln_sss         !: Logical switch for sea surface salinity 
     84   LOGICAL, PUBLIC :: ln_sstnight    !: Logical switch for night mean SST observations 
    8285   LOGICAL, PUBLIC :: ln_nea         !: Remove observations near land 
    8386   LOGICAL, PUBLIC :: ln_altbias     !: Logical switch for altimeter bias   
    … …  
    167170         &            nmsshc, mdtcorr, mdtcutoff,                     & 
    168171         &            ln_reysst, ln_ghrsst, reysstname, reysstfmt,    & 
     172         &            ln_sstnight,                                    & 
    169173         &            ln_grid_search_lookup,                          & 
    170174         &            grid_search_file, grid_search_res,              & 
    … …  
    176180         &            ln_velhradcp, velhradcpfiles,                   & 
    177181         &            ln_velfb, velfbfiles, ln_velfb_av,              & 
    178          &            ln_profb_enatim, ln_ignmis 
     182         &            ln_profb_enatim, ln_ignmis, ln_cl4 
    179183 
    180184      INTEGER :: jprofset 
    … …  
    226230      ln_velhradcp = .FALSE. 
    227231      ln_velfb = .FALSE. 
     232      ln_sstnight = .FALSE. 
    228233      ln_nea = .FALSE. 
    229234      ln_grid_search_lookup = .FALSE. 
    230235      ln_grid_global = .FALSE. 
    231236      ln_s_at_t = .TRUE. 
     237      ln_cl4 = .FALSE. 
    232238      grid_search_file = 'xypos' 
    233239      bias_file='bias.nc' 
    … …  
    357363         WRITE(numout,*) '             Logical switch for GHRSST observations          ln_ghrsst = ', ln_ghrsst 
    358364         WRITE(numout,*) '             Logical switch for feedback SST data             ln_sstfb = ', ln_sstfb 
     365         WRITE(numout,*) '             Logical switch for night-time SST obs         ln_sstnight = ', ln_sstnight 
    359366         WRITE(numout,*) '             Logical switch for SSS observations                ln_sss = ', ln_sss 
    360367         WRITE(numout,*) '             Logical switch for Sea Ice observations         ln_seaice = ', ln_seaice 
    … …  
    750757            nsstsets = nsstsets + 1 
    751758 
    752             ld_sstnight(nsstsets) = .TRUE. 
     759            ld_sstnight(nsstsets) = ln_sstnight 
    753760 
    754761            CALL obs_rea_sst_rey( reysstname, reysstfmt, sstdata(nsstsets), & 
    … …  
    764771            nsstsets = nsstsets + 1 
    765772 
    766             ld_sstnight(nsstsets) = .TRUE. 
     773            ld_sstnight(nsstsets) = ln_sstnight 
    767774           
    768775            CALL obs_rea_sst( 1, sstdata(nsstsets), jnumsst, & 
    … …  
    783790               nsstsets = nsstsets + 1 
    784791 
    785                ld_sstnight(nsstsets) = .TRUE. 
     792               ld_sstnight(nsstsets) = ln_sstnight 
    786793             
    787794               CALL obs_rea_sst( 0, sstdata(nsstsets), 1, & 
    … …  
    14141421   END SUBROUTINE dia_obs_wri 
    14151422 
     1423   SUBROUTINE dia_obs_dealloc 
     1424      IMPLICIT NONE 
     1425      !!---------------------------------------------------------------------- 
     1426      !!                    *** ROUTINE dia_obs_dealloc *** 
     1427      !! 
     1428      !!  ** Purpose : To deallocate data to enable the obs_oper online loop. 
     1429      !!               Specifically: dia_obs_init --> dia_obs --> dia_obs_wri 
     1430      !! 
     1431      !!  ** Method : Clean up various arrays left behind by the obs_oper. 
     1432      !! 
     1433      !!  ** Action : 
     1434      !! 
     1435      !!---------------------------------------------------------------------- 
     1436      !! obs_grid deallocation 
     1437      CALL obs_grid_deallocate 
     1438 
     1439      !! diaobs deallocation 
     1440      IF ( nprofsets > 0 ) THEN 
     1441          DEALLOCATE(ld_enact, & 
     1442                  &  profdata, & 
     1443                  &  prodatqc) 
     1444      END IF 
     1445      IF ( ln_sla ) THEN 
     1446          DEALLOCATE(sladata, & 
     1447                  &  sladatqc) 
     1448      END IF 
     1449      IF ( ln_seaice ) THEN 
     1450          DEALLOCATE(sladata, & 
     1451                  &  sladatqc) 
     1452      END IF 
     1453      IF ( ln_sst ) THEN 
     1454          DEALLOCATE(sstdata, & 
     1455                  &  sstdatqc) 
     1456      END IF 
     1457      IF ( ln_vel3d ) THEN 
     1458          DEALLOCATE(ld_velav, & 
     1459                  &  velodata, & 
     1460                  &  veldatqc) 
     1461      END IF 
     1462   END SUBROUTINE dia_obs_dealloc 
     1463 
    14161464   SUBROUTINE ini_date( ddobsini ) 
    14171465      !!---------------------------------------------------------------------- 

  • branches/2013/dev_UKMO_2013/NEMOGCM/NEMO/OPA_SRC/OBS/mpp_map.F90

    r2363 r4171  
    5252      !!---------------------------------------------------------------------- 
    5353 
    54       ALLOCATE( & 
    55          & mppmap(jpiglo,jpjglo) & 
    56          & ) 
    57  
     54      IF (.NOT. ALLOCATED(mppmap)) THEN 
     55         ALLOCATE( & 
     56            & mppmap(jpiglo,jpjglo) & 
     57            & ) 
     58      ENDIF 
    5859      ! Initialize local imppmap 
    5960 

  • branches/2013/dev_UKMO_2013/NEMOGCM/NEMO/OPA_SRC/OBS/obs_fbm.F90

    r2287 r4171  
    4545   INTEGER, PARAMETER    :: fbimdi = -99999   !: Integers 
    4646   REAL(fbsp), PARAMETER :: fbrmdi =  99999   !: Reals 
    47     
     47 
     48   ! Output stream choice 
     49   LOGICAL               :: ln_cl4 = .FALSE.  !: Logical switch for 
     50                                              !: class 4 file outputs 
     51  
    4852   ! Main data structure for observation feedback data. 
    4953 
    … …  
    10261030 
    10271031   SUBROUTINE write_obfbdata( cdfilename, fbdata ) 
     1032      !!---------------------------------------------------------------------- 
     1033      !!                    ***  ROUTINE write_obfbdata  *** 
     1034      !! 
     1035      !! ** Purpose :   Write an obfbdata structure into a netCDF file. 
     1036      !! 
     1037      !! ** Method  :   Decides which output wrapper to use.  
     1038      !! 
     1039      !! ** Action  :  
     1040      !! 
     1041      !!---------------------------------------------------------------------- 
     1042      !! * Arguments 
     1043      CHARACTER(len=*) :: cdfilename ! Output filename 
     1044      TYPE(obfbdata)   :: fbdata     ! obsfbdata structure 
     1045#if defined key_offobsoper 
     1046      IF (ln_cl4) THEN 
     1047          ! Class 4 file output stream 
     1048          CALL write_obfbdata_cl( cdfilename, fbdata ) 
     1049      ELSE 
     1050#endif 
     1051          ! Standard feedback file output stream 
     1052          CALL write_obfbdata_fb( cdfilename, fbdata ) 
     1053#if defined key_offobsoper 
     1054      ENDIF 
     1055#endif 
     1056   END SUBROUTINE write_obfbdata 
     1057 
     1058   SUBROUTINE write_obfbdata_fb( cdfilename, fbdata ) 
    10281059      !!---------------------------------------------------------------------- 
    10291060      !!                    ***  ROUTINE write_obfbdata  *** 
    … …  
    15241555 
    15251556       
    1526    END SUBROUTINE write_obfbdata 
     1557   END SUBROUTINE write_obfbdata_fb 
     1558 
     1559#if defined key_offobsoper 
     1560   SUBROUTINE write_obfbdata_cl(cdfilename, fbdata) 
     1561      !!---------------------------------------------------------------------- 
     1562      !!                    ***  ROUTINE write_obfbdata_cl  *** 
     1563      !! 
     1564      !! ** Purpose : Write an obfbdata structure into a class 4 file. 
     1565      !! 
     1566      !! ** Method  : 1. Allocate memory needed by ooo_write 
     1567      !!              2. Map obfbdata into allocated memory 
     1568      !!              3. Pass mapped data to ooo_write 
     1569      !!              4. Deallocate memory 
     1570      !!---------------------------------------------------------------------- 
     1571      USE dom_oce, ONLY: narea 
     1572      USE ooo_write 
     1573      USE ooo_data 
     1574      !! * Arguments 
     1575      CHARACTER(len=*) :: cdfilename ! Feedback filename 
     1576      TYPE(obfbdata)   :: fbdata     ! obsfbdata structure 
     1577      !! * Local variables 
     1578      CHARACTER(len=17), PARAMETER :: cpname = 'write_obfbdata_cl' 
     1579      CHARACTER(len=64) :: & 
     1580              & cdate, &   !: class 4 file validity date  
     1581              & cconf, &   !: model configuration 
     1582              & csys, &    !: model system 
     1583              & ccont, &   !: contact email 
     1584              & cinst, &   !: institution 
     1585              & cversion   !: model version 
     1586      CHARACTER(len=8) :: & 
     1587              & ckind      !: observation kind 
     1588      CHARACTER(len=3) :: cfield 
     1589      INTEGER :: kobs, &   !: number of observations 
     1590              &  kvars, &  !: number of physical variables 
     1591              &  kdeps, &  !: number of observed depths 
     1592              &  kfcst, &  !: number of forecasts 
     1593              &  kifcst, & !: current forecast number 
     1594              &  kproc     !: processor number 
     1595      INTEGER, DIMENSION(:, :, :), ALLOCATABLE :: & 
     1596              &  kqc       !: quality control counterpart 
     1597      INTEGER(KIND=2), DIMENSION(:, :, :), ALLOCATABLE :: & 
     1598              &  k2qc       !: quality control counterpart 
     1599      REAL(kind=fbdp) :: & 
     1600              &  pmodjuld  !: model Julian day 
     1601      REAL(kind=fbdp), DIMENSION(:), ALLOCATABLE :: & 
     1602              &  plead, &  !: forecast lead time 
     1603              &  plam, &   !: longitude of observation 
     1604              &  pphi, &   !: latitude of observation 
     1605              &  ptim      !: time of observation 
     1606      REAL(kind=fbdp), DIMENSION(:, :), ALLOCATABLE :: & 
     1607              &  pdep      !: depths of observations 
     1608      REAL(kind=fbdp), DIMENSION(:, :, :), ALLOCATABLE :: & 
     1609              &  pob, &    !: observation counterpart 
     1610              &  pextra    !: extra field counterpart 
     1611      REAL(kind=fbdp), DIMENSION(:, :, :), ALLOCATABLE :: & 
     1612              &  pmod      !: model counterpart 
     1613      CHARACTER(len=128) :: & 
     1614              &  clfilename  !: class 4 file name 
     1615      CHARACTER(len=128), DIMENSION(:), ALLOCATABLE :: & 
     1616              &  ctype       !: Instrument type 
     1617      CHARACTER(len=nf90_max_name) :: & 
     1618              & cdtmp        !: NetCDF variable name 
     1619      CHARACTER(len=8), DIMENSION(:), ALLOCATABLE :: & 
     1620              &  cwmo, &     !: Instrument WMO ID 
     1621              &  cunit, &    !: Instrument WMO ID 
     1622              &  cvarname    !: Instrument WMO ID 
     1623      INTEGER :: & 
     1624              &  idep, &     !: Loop variable 
     1625              &  ivar, &     !: Loop variable 
     1626              &  iobs, &     !: Loop variable 
     1627              &  ii, &       !: Loop variable 
     1628              &  ij, &       !: Loop variable 
     1629              &  ik, &       !: Loop variable 
     1630              &  il          !: Loop variable 
     1631      cconf = TRIM(cl4_cfg) 
     1632      csys = TRIM(cl4_sys) 
     1633      cversion = TRIM(cl4_vn) 
     1634      ccont = TRIM(cl4_contact) 
     1635      cinst = TRIM(cl4_inst) 
     1636      cdate = TRIM(cl4_date) 
     1637      CALL locate_kind(cdfilename, ckind) 
     1638      kproc = narea 
     1639      kfcst = cl4_fcst_len 
     1640      kobs = fbdata%nobs 
     1641      kdeps = fbdata%nlev 
     1642      kvars = fbdata%nvar 
     1643      IF (kobs .GT. 0) THEN 
     1644         ALLOCATE(plam(kobs), & 
     1645               &  pphi(kobs), & 
     1646               &  ptim(kobs), & 
     1647               &  plead(kfcst), & 
     1648               &  pdep(kdeps, kobs), & 
     1649               &  kqc(kdeps, kvars, kobs), & 
     1650               &  k2qc(kdeps, kvars, kobs), & 
     1651               &  pob(kdeps, kvars, kobs), & 
     1652               &  pmod(kdeps, kvars, kobs), & 
     1653               &  pextra(kdeps, kvars, kobs), & 
     1654               &  ctype(kobs), & 
     1655               &  cwmo(kobs), & 
     1656               &  cunit(kvars), & 
     1657               &  cvarname(kvars)) 
     1658         plam(:) = fbdata%plam(:) 
     1659         pphi(:) = fbdata%pphi(:) 
     1660         ptim(:) = fbdata%ptim(:) 
     1661         pdep(:, :) = fbdata%pdep(:, :) 
     1662         kqc(:,:,:) = 1. 
     1663         DO ii = 1, kvars 
     1664            cvarname(ii)  = fbdata%cname(ii) 
     1665            cunit(ii)     = fbdata%cobunit(ii) 
     1666         END DO 
     1667 
     1668         ! Quality control algorithm 
     1669         k2qc(:,:,:) = NF90_FILL_SHORT 
     1670         DO idep = 1,kdeps 
     1671            DO ivar = 1, kvars 
     1672               DO iobs = 1, kobs 
     1673                  ! 1 symbolises good for fbdata 
     1674                  ! fbimdi symbolises that qc has not been set 
     1675                  ! Essentially, if any fbdata flag is not an element of {1, fbimdi} 
     1676                  ! then set the class 4 flag to bad. 
     1677                  ! Note: fbdata%ioqc is marked good if zero. 
     1678                  IF (((fbdata%ioqc(iobs) /= 0) .AND. & 
     1679                            & (fbdata%ioqc(iobs) /= fbimdi)) .OR. & 
     1680                    & ((fbdata%ipqc(iobs) /= 1) .AND. & 
     1681                            & (fbdata%ipqc(iobs) /= fbimdi)) .OR. & 
     1682                    & ((fbdata%idqc(idep,iobs) /= 1) .AND. & 
     1683                            & (fbdata%idqc(idep,iobs) /= fbimdi)) .OR. & 
     1684                    & ((fbdata%ivqc(iobs,ivar) /= 1) .AND. & 
     1685                            & (fbdata%ivqc(iobs,ivar) /= fbimdi)) .OR. & 
     1686                    & ((fbdata%ivlqc(idep,iobs,ivar) /= 1) .AND. & 
     1687                            & (fbdata%ivlqc(idep,iobs,ivar) /= fbimdi)) .OR. & 
     1688                    & ((fbdata%itqc(iobs) /= 1) .AND. & 
     1689                            & (fbdata%itqc(iobs) /= fbimdi))) THEN 
     1690                     ! 1 symbolises bad for class 4 file 
     1691                     k2qc(idep, ivar, iobs) = 1 
     1692                  ELSE 
     1693                     ! 0 symbolises good for class 4 file 
     1694                     k2qc(idep, ivar, iobs) = 0 
     1695                  END IF  
     1696               END DO 
     1697            END DO 
     1698         END DO 
     1699 
     1700         ! Permute observation dimensions 
     1701         pob(:,:,:) = RESHAPE(fbdata%pob, (/kdeps, kvars, kobs/), & 
     1702                            & ORDER=(/1, 3, 2/)) 
     1703 
     1704         ! Explicit model counterpart dimension permutation 
     1705         ! 1,2,3,4 --> 1,4,2,3 
     1706         pmod(:,:,:) = fbrmdi 
     1707         ij = cl4_fcst_idx(jimatch) 
     1708         DO ii = 1,kdeps 
     1709            DO ik = 1, kvars 
     1710               DO il = 1, kobs 
     1711                  pmod(ii,ik,il) = fbdata%padd(ii,il,1,ik) 
     1712               END DO 
     1713            END DO 
     1714         END DO 
     1715 
     1716         ! Extra fields set to missing for now 
     1717         pextra(:,:,:) = fbrmdi 
     1718 
     1719         ! Lead time of class 4 file is a global parameter 
     1720         plead = cl4_leadtime(1:cl4_fcst_len) 
     1721 
     1722         ! Model Julian day 
     1723         pmodjuld = cl4_modjuld 
     1724 
     1725         ! Observation types 
     1726         ctype(:) = 'X' 
     1727         DO ii = 1,kobs 
     1728            ctype(ii) = fbdata%cdtyp(ii) 
     1729         END DO 
     1730 
     1731         ! World Meteorology Organisation codes 
     1732         cwmo(:) = fbdata%cdwmo(:) 
     1733 
     1734         ! Initialise class 4 file 
     1735         CALL ooo_wri_init(cconf, csys, ckind, cversion, ccont, cinst, cdate, & 
     1736                         & kproc, kobs, kvars, kdeps, kfcst, & 
     1737                         & clfilename) 
     1738 
     1739         ! Write standard variables 
     1740         CALL ooo_wri_default(clfilename, kobs, kvars, kfcst, kdeps, & 
     1741                            & ctype, cwmo, cunit, cvarname, & 
     1742                            & plam, pphi, pdep, ptim, pob, plead, & 
     1743                            & k2qc, pmodjuld) 
     1744         !! Write to optional variables 
     1745         cdtmp = cl4_vars(jimatch) 
     1746         IF ( (TRIM(cdtmp) == "forecast") .OR. & 
     1747              (TRIM(cdtmp) == "persistence") ) THEN 
     1748            !! 4D variables 
     1749            CALL ooo_wri_extra(clfilename, TRIM(cdtmp), kdeps, kfcst, & 
     1750                            &  kvars, kobs, (/ 1,ij,1,1 /), (/ kdeps,1,kvars,kobs /), pmod) 
     1751         ELSE 
     1752            !! 3D variables 
     1753            CALL ooo_wri_extra(clfilename, TRIM(cdtmp), kdeps, & 
     1754                            &  kvars, kobs, (/ 1,1,1 /), (/ kdeps,kvars,kobs /), pmod) 
     1755         ENDIF 
     1756 
     1757         DEALLOCATE(plam, pphi, ptim, pdep, plead, kqc, k2qc, & 
     1758                  & pob, pmod, pextra, ctype, cwmo, & 
     1759                  & cunit, cvarname) 
     1760      END IF 
     1761   END SUBROUTINE write_obfbdata_cl 
     1762#endif 
     1763 
     1764#if defined key_offobsoper 
     1765   SUBROUTINE locate_kind(cdfilename, ckind) 
     1766      !!---------------------------------------------------------------------- 
     1767      !!                    ***  ROUTINE locate_kind  *** 
     1768      !! 
     1769      !! ** Purpose : Detect which kind of class 4 file is being produced. 
     1770      !! 
     1771      !! ** Method  : 1. Inspect cdfilename for observation kind. 
     1772      !!---------------------------------------------------------------------- 
     1773      CHARACTER(len=*) :: cdfilename ! Feedback filename 
     1774      CHARACTER(len=8) :: ckind 
     1775      IF (cdfilename(1:3) == 'sst') THEN 
     1776         ckind = 'SST' 
     1777      ELSE IF (cdfilename(1:3) == 'sla') THEN 
     1778         ckind = 'SLA' 
     1779      ELSE IF (cdfilename(1:3) == 'pro') THEN 
     1780         ckind = 'profile' 
     1781      ELSE IF (cdfilename(1:3) == 'ena') THEN 
     1782         ckind = 'profile' 
     1783      ELSE IF (cdfilename(1:3) == 'sea') THEN 
     1784         ckind = 'seaice' 
     1785      ELSE 
     1786         ckind = 'unknown' 
     1787      END IF 
     1788   END SUBROUTINE locate_kind 
     1789#endif 
    15271790 
    15281791   SUBROUTINE putvaratt_obfbdata( idfile, idvar, cdlongname, cdunits, & 

  • branches/2013/dev_UKMO_2013/NEMOGCM/NEMO/OPA_SRC/OBS/obs_oper.F90

    r3651 r4171  
    861861 
    862862         ENDIF 
    863           
    864863         sstdatqc%rmod(jobs,1) = zext(1) 
    865864          

  • branches/2013/dev_UKMO_2013/NEMOGCM/NEMO/OPA_SRC/SBC/geo2ocean.F90

    r2715 r4171  
    187187         &      gsinf(jpi,jpj), gcosf(jpi,jpj), STAT=ierr ) 
    188188      IF(lk_mpp)   CALL mpp_sum( ierr ) 
    189       IF( ierr /= 0 )   CALL ctl_stop('STOP', 'angle_msh_geo: unable to allocate arrays' ) 
     189      IF( ierr /= 0 )   CALL ctl_stop('angle: unable to allocate arrays' ) 
    190190 
    191191      ! ============================= ! 
    … …  
    361361            &      gsinlat(jpi,jpj,4) , gcoslat(jpi,jpj,4) , STAT=ierr ) 
    362362         IF( lk_mpp    )   CALL mpp_sum( ierr ) 
    363          IF( ierr /= 0 )   CALL ctl_stop('STOP', 'angle_msh_geo: unable to allocate arrays' ) 
     363         IF( ierr /= 0 )   CALL ctl_stop('geo2oce: unable to allocate arrays' ) 
    364364      ENDIF 
    365365 
    … …  
    438438      !!---------------------------------------------------------------------- 
    439439 
    440       IF( ALLOCATED( gsinlon ) ) THEN 
     440      IF( .NOT. ALLOCATED( gsinlon ) ) THEN 
    441441         ALLOCATE( gsinlon(jpi,jpj,4) , gcoslon(jpi,jpj,4) ,   & 
    442442            &      gsinlat(jpi,jpj,4) , gcoslat(jpi,jpj,4) , STAT=ierr ) 
    443443         IF( lk_mpp    )   CALL mpp_sum( ierr ) 
    444          IF( ierr /= 0 )   CALL ctl_stop('STOP', 'angle_msh_geo: unable to allocate arrays' ) 
     444         IF( ierr /= 0 )   CALL ctl_stop('oce2geo: unable to allocate arrays' ) 
    445445      ENDIF 
    446446 

  • branches/2013/dev_UKMO_2013/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90

    r3914 r4171  
    388388      ! 
    389389      IF( TRIM( sn_rcv_tau%cldes ) /= 'oce and ice' ) THEN        ! 'oce and ice' case ocean stress on ocean mesh used 
    390          srcv(jpr_itz1:jpr_itz2)%laction = .FALSE.    ! ice components not received (itx1 and ity1 used later) 
     390         srcv(jpr_itx1:jpr_itz2)%laction = .FALSE.    ! ice components not received 
    391391         srcv(jpr_itx1)%clgrid = 'U'                  ! ocean stress used after its transformation 
    392392         srcv(jpr_ity1)%clgrid = 'V'                  ! i.e. it is always at U- & V-points for i- & j-comp. resp. 
    … …  
    407407      SELECT CASE( TRIM( sn_rcv_emp%cldes ) ) 
    408408      CASE( 'oce only'      )   ;   srcv(                                 jpr_oemp   )%laction = .TRUE.  
    409       CASE( 'conservative'  )   ;   srcv( (/jpr_rain, jpr_snow, jpr_ievp, jpr_tevp/) )%laction = .TRUE. 
     409      CASE( 'conservative'  ) 
     410         srcv( (/jpr_rain, jpr_snow, jpr_ievp, jpr_tevp/) )%laction = .TRUE. 
     411         IF ( k_ice <= 1 )  srcv(jpr_ivep)%laction = .FALSE. 
    410412      CASE( 'oce and ice'   )   ;   srcv( (/jpr_ievp, jpr_sbpr, jpr_semp, jpr_oemp/) )%laction = .TRUE. 
    411413      CASE default              ;   CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_rcv_emp%cldes' ) 
    … …  
    465467         CALL ctl_stop( 'sbc_cpl_init: namsbc_cpl namelist mismatch between sn_rcv_qns%cldes and sn_rcv_dqnsdt%cldes' ) 
    466468      !                                                      ! ------------------------- ! 
    467       !                                                      !    Ice Qsr penetration    !    
    468       !                                                      ! ------------------------- ! 
    469       ! fraction of net shortwave radiation which is not absorbed in the thin surface layer  
    470       ! and penetrates inside the ice cover ( Maykut and Untersteiner, 1971 ; Elbert anbd Curry, 1993 ) 
    471       ! Coupled case: since cloud cover is not received from atmosphere  
    472       !               ===> defined as constant value -> definition done in sbc_cpl_init 
    473       fr1_i0(:,:) = 0.18 
    474       fr2_i0(:,:) = 0.82 
    475       !                                                      ! ------------------------- ! 
    476469      !                                                      !      10m wind module      !    
    477470      !                                                      ! ------------------------- ! 
    … …  
    508501      ! Allocate taum part of frcv which is used even when not received as coupling field 
    509502      IF ( .NOT. srcv(jpr_taum)%laction ) ALLOCATE( frcv(jpr_taum)%z3(jpi,jpj,srcv(jn)%nct) ) 
     503      ! Allocate itx1 and ity1 as they are used in sbc_cpl_ice_tau even if srcv(jpr_itx1)%laction = .FALSE. 
     504      IF( k_ice /= 0 ) THEN 
     505         IF ( .NOT. srcv(jpr_itx1)%laction ) ALLOCATE( frcv(jpr_itx1)%z3(jpi,jpj,srcv(jn)%nct) ) 
     506         IF ( .NOT. srcv(jpr_ity1)%laction ) ALLOCATE( frcv(jpr_ity1)%z3(jpi,jpj,srcv(jn)%nct) ) 
     507      END IF 
    510508 
    511509      ! ================================ ! 
    … …  
    13291327      END SELECT 
    13301328 
     1329      !    Ice Qsr penetration used (only?)in lim2 or lim3  
     1330      ! fraction of net shortwave radiation which is not absorbed in the thin surface layer  
     1331      ! and penetrates inside the ice cover ( Maykut and Untersteiner, 1971 ; Elbert anbd Curry, 1993 ) 
     1332      ! Coupled case: since cloud cover is not received from atmosphere  
     1333      !               ===> defined as constant value -> definition done in sbc_cpl_init 
     1334      fr1_i0(:,:) = 0.18 
     1335      fr2_i0(:,:) = 0.82 
     1336 
     1337 
    13311338      CALL wrk_dealloc( jpi,jpj, zcptn, ztmp, zicefr ) 
    13321339      ! 

  • branches/2013/dev_UKMO_2013/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90

    r3905 r4171  
    221221      ENDIF 
    222222      ! 
     223                          CALL sbc_ssm_init               ! Sea-surface mean fields initialisation 
     224      ! 
    223225      IF( ln_ssr      )   CALL sbc_ssr_init               ! Sea-Surface Restoring initialisation 
    224226      ! 

  • branches/2013/dev_UKMO_2013/NEMOGCM/NEMO/OPA_SRC/TRA/eosbn2.F90

    r3625 r4171  
    675675 
    676676 
    677    FUNCTION tfreez( psal ) RESULT( ptf ) 
     677   FUNCTION tfreez( psal, pdep ) RESULT( ptf ) 
    678678      !!---------------------------------------------------------------------- 
    679679      !!                 ***  ROUTINE eos_init  *** 
    … …  
    688688      !!---------------------------------------------------------------------- 
    689689      REAL(wp), DIMENSION(jpi,jpj), INTENT(in   ) ::   psal   ! salinity             [psu] 
     690      REAL(wp), DIMENSION(jpi,jpj), INTENT(in   ), OPTIONAL ::   pdep   ! depth      [decibars] 
    690691      ! Leave result array automatic rather than making explicitly allocated 
    691692      REAL(wp), DIMENSION(jpi,jpj)                ::   ptf    ! freezing temperature [Celcius] 
    … …  
    694695      ptf(:,:) = ( - 0.0575_wp + 1.710523e-3_wp * SQRT( psal(:,:) )   & 
    695696         &                     - 2.154996e-4_wp *       psal(:,:)   ) * psal(:,:) 
     697      IF ( PRESENT( pdep ) ) THEN    
     698         ptf(:,:) = ptf(:,:) - 7.53e-4_wp * pdep(:,:) 
     699      ENDIF 
    696700      ! 
    697701   END FUNCTION tfreez 

  • branches/2013/dev_UKMO_2013/NEMOGCM/NEMO/SAS_SRC/daymod.F90

    r3851 r4171  
    246246               nday_year = 1 
    247247               nsec_year = ndt05 
     248               IF( nsec1jan000 >= 2 * (2**30 - nsecd * nyear_len(1) / 2 ) ) THEN   ! test integer 4 max value 
     249                  CALL ctl_stop( 'The number of seconds between Jan. 1st 00h of nit000 year and Jan. 1st 00h ',   & 
     250                     &           'of the current year is exceeding the INTEGER 4 max VALUE: 2^31-1 -> 68.09 years in seconds', & 
     251                     & 'You must do a restart at higher frequency (or remove this STOP and recompile everything in I8)' ) 
     252               ENDIF 
    248253               nsec1jan000 = nsec1jan000 + nsecd * nyear_len(1) 
    249254               IF( nleapy == 1 )   CALL day_mth 

  • branches/2013/dev_UKMO_2013/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zsed.F90

    r3905 r4171  
    8282      IF( nn_timing == 1 )  CALL timing_start('p4z_sed') 
    8383      ! 
    84       IF( kt == nit000 .AND. jnt == 1 )  THEN 
     84      IF( kt == nittrc000 .AND. jnt == 1 )  THEN 
    8585         ryyss    = nyear_len(1) * rday    ! number of seconds per year and per month 
    8686         rmtss    = ryyss / raamo 

  • branches/2013/dev_UKMO_2013/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zsms.F90

    r3972 r4171  
    7676      ENDIF 
    7777      ! 
    78       IF( ln_rsttr .AND. kt == nittrc000 )                         CALL p4z_rst( nittrc000, 'READ' )  !* read or initialize all required fields  
     78      IF( kt == nittrc000 ) THEN 
     79        ! 
     80        CALL p4z_che                              ! initialize the chemical constants 
     81        ! 
     82        IF( .NOT. ln_rsttr ) THEN  ;   CALL p4z_ph_ini   !  set PH at kt=nit000  
     83        ELSE                       ;   CALL p4z_rst( nittrc000, 'READ' )  !* read or initialize all required fields  
     84        ENDIF 
     85        ! 
     86      ENDIF 
     87 
    7988      IF( ln_pisdmp .AND. MOD( kt - nn_dttrc, nn_pisdmp ) == 0 )   CALL p4z_dmp( kt )      ! Relaxation of some tracers 
    8089      ! 
    … …  
    238247   END SUBROUTINE p4z_sms_init 
    239248 
     249   SUBROUTINE p4z_ph_ini 
     250      !!--------------------------------------------------------------------- 
     251      !!                   ***  ROUTINE p4z_ini_ph  *** 
     252      !! 
     253      !!  ** Purpose : Initialization of chemical variables of the carbon cycle 
     254      !!--------------------------------------------------------------------- 
     255      INTEGER  ::  ji, jj, jk 
     256      REAL(wp) ::  zcaralk, zbicarb, zco3 
     257      REAL(wp) ::  ztmas, ztmas1 
     258      !!--------------------------------------------------------------------- 
     259 
     260      ! Set PH from  total alkalinity, borat (???), akb3 (???) and ak23 (???) 
     261      ! -------------------------------------------------------- 
     262      DO jk = 1, jpk 
     263         DO jj = 1, jpj 
     264            DO ji = 1, jpi 
     265               ztmas   = tmask(ji,jj,jk) 
     266               ztmas1  = 1. - tmask(ji,jj,jk) 
     267               zcaralk = trn(ji,jj,jk,jptal) - borat(ji,jj,jk) / (  1. + 1.E-8 / ( rtrn + akb3(ji,jj,jk) )  ) 
     268               zco3    = ( zcaralk - trn(ji,jj,jk,jpdic) ) * ztmas + 0.5e-3 * ztmas1 
     269               zbicarb = ( 2. * trn(ji,jj,jk,jpdic) - zcaralk ) 
     270               hi(ji,jj,jk) = ( ak23(ji,jj,jk) * zbicarb / zco3 ) * ztmas + 1.e-9 * ztmas1 
     271            END DO 
     272         END DO 
     273     END DO 
     274     ! 
     275   END SUBROUTINE p4z_ph_ini 
     276 
    240277   SUBROUTINE p4z_rst( kt, cdrw ) 
    241278      !!--------------------------------------------------------------------- 
    … …  
    266303         ELSE 
    267304!            hi(:,:,:) = 1.e-9  
    268             ! Set PH from  total alkalinity, borat (???), akb3 (???) and ak23 (???) 
    269             ! -------------------------------------------------------- 
    270             DO jk = 1, jpk 
    271                DO jj = 1, jpj 
    272                   DO ji = 1, jpi 
    273                      ztmas   = tmask(ji,jj,jk) 
    274                      ztmas1  = 1. - tmask(ji,jj,jk) 
    275                      zcaralk = trn(ji,jj,jk,jptal) - borat(ji,jj,jk) / (  1. + 1.E-8 / ( rtrn + akb3(ji,jj,jk) )  ) 
    276                      zco3    = ( zcaralk - trn(ji,jj,jk,jpdic) ) * ztmas + 0.5e-3 * ztmas1 
    277                      zbicarb = ( 2. * trn(ji,jj,jk,jpdic) - zcaralk ) 
    278                      hi(ji,jj,jk) = ( ak23(ji,jj,jk) * zbicarb / zco3 ) * ztmas + 1.e-9 * ztmas1 
    279                   END DO 
    280                END DO 
    281             END DO 
     305            CALL p4z_ph_ini 
    282306         ENDIF 
    283307         CALL iom_get( numrtr, jpdom_autoglo, 'Silicalim', xksi(:,:) ) 

  • branches/2013/dev_UKMO_2013/NEMOGCM/NEMO/TOP_SRC/PISCES/trcini_pisces.F90

    r3757 r4171  
    122122      rdenita =   3._wp /  5._wp 
    123123      o2ut    = 131._wp / 122._wp 
    124  
    125       CALL p4z_che        ! initialize the chemical constants 
    126124 
    127125      ! Initialization of tracer concentration in case of  no restart  
    … …  
    162160         xksi(:,:)    = 2.e-6 
    163161         xksimax(:,:) = xksi(:,:) 
    164  
    165          ! Initialization of chemical variables of the carbon cycle 
    166          ! -------------------------------------------------------- 
    167          DO jk = 1, jpk 
    168             DO jj = 1, jpj 
    169                DO ji = 1, jpi 
    170                   ztmas   = tmask(ji,jj,jk) 
    171                   ztmas1  = 1. - tmask(ji,jj,jk) 
    172                   zcaralk = trn(ji,jj,jk,jptal) - borat(ji,jj,jk) / (  1. + 1.E-8 / ( rtrn + akb3(ji,jj,jk) )  ) 
    173                   zco3    = ( zcaralk - trn(ji,jj,jk,jpdic) ) * ztmas + 0.5e-3 * ztmas1 
    174                   zbicarb = ( 2. * trn(ji,jj,jk,jpdic) - zcaralk ) 
    175                   hi(ji,jj,jk) = ( ak23(ji,jj,jk) * zbicarb / zco3 ) * ztmas + 1.e-9 * ztmas1 
    176                END DO 
    177             END DO 
    178          END DO 
    179          ! 
     162        ! 
    180163      END IF 
    181164 

  • branches/2013/dev_UKMO_2013/NEMOGCM/TOOLS/COMPILE/Fcheck_archfile.sh

    r3925 r4171  
    4040# :: 
    4141# 
    42 #  $ ./Fcheck_archfile.sh ARCHFILE COMPILER 
     42#  $ ./Fcheck_archfile.sh ARCHFILE CPPFILE COMPILER 
    4343# 
    4444# 
    … …  
    9494   else 
    9595       if [ -f ${COMPIL_DIR}/$1 ]; then 
    96       # has the cpp keys file been changed since we copied the arch file in ${COMPIL_DIR}? 
    97       mycpp=$( ls -l ${COMPIL_DIR}/$2 | sed -e "s/.* -> //" ) 
    98       if [ "$mycpp" != "$( cat ${COMPIL_DIR}/cpp.history )" ]; then 
    99           echo $mycpp > ${COMPIL_DIR}/cpp.history 
    100           cpeval ${myarch} ${COMPIL_DIR}/$1 
     96      if [ "$2" != "nocpp" ]  
     97      then 
     98          # has the cpp keys file been changed since we copied the arch file in ${COMPIL_DIR}? 
     99          mycpp=$( ls -l ${COMPIL_DIR}/$2 | sed -e "s/.* -> //" ) 
     100          if [ "$mycpp" != "$( cat ${COMPIL_DIR}/cpp.history )" ]; then 
     101         echo $mycpp > ${COMPIL_DIR}/cpp.history 
     102         cpeval ${myarch} ${COMPIL_DIR}/$1 
     103          fi 
     104          # has the cpp keys file been updated since we copied the arch file in ${COMPIL_DIR}? 
     105          mycpp=$( find -L ${COMPIL_DIR} -cnewer ${COMPIL_DIR}/$1 -name $2 -print ) 
     106          [ ${#mycpp} -ne 0 ] && cpeval ${myarch} ${COMPIL_DIR}/$1 
    101107      fi 
    102       # has the cpp keys file been updated since we copied the arch file in ${COMPIL_DIR}? 
    103       mycpp=$( find -L ${COMPIL_DIR} -cnewer ${COMPIL_DIR}/$1 -name $2 -print ) 
    104       [ ${#mycpp} -ne 0 ] && cpeval ${myarch} ${COMPIL_DIR}/$1 
    105108      # has myarch file been updated since we copied it in ${COMPIL_DIR}? 
    106109      myarchdir=$( dirname ${myarch} ) 
    … …  
    134137    if [ "$myarch" == "$( cat ${COMPIL_DIR}/arch.history )" ]; then  
    135138   if [ -f ${COMPIL_DIR}/$1 ]; then 
    136        # has the cpp keys file been changed since we copied the arch file in ${COMPIL_DIR}? 
    137        mycpp=$( ls -l ${COMPIL_DIR}/$2 | sed -e "s/.* -> //" ) 
    138        if [ "$mycpp" != "$( cat ${COMPIL_DIR}/cpp.history )" ]; then 
    139       echo $mycpp > ${COMPIL_DIR}/cpp.history 
    140       cpeval ${myarch} ${COMPIL_DIR}/$1 
     139       if [ "$2" != "nocpp" ]  
     140       then 
     141      # has the cpp keys file been changed since we copied the arch file in ${COMPIL_DIR}? 
     142      mycpp=$( ls -l ${COMPIL_DIR}/$2 | sed -e "s/.* -> //" ) 
     143      if [ "$mycpp" != "$( cat ${COMPIL_DIR}/cpp.history )" ]; then 
     144          echo $mycpp > ${COMPIL_DIR}/cpp.history 
     145          cpeval ${myarch} ${COMPIL_DIR}/$1 
     146      fi 
     147      # has the cpp keys file been updated since we copied the arch file in ${COMPIL_DIR}? 
     148      mycpp=$( find -L ${COMPIL_DIR} -cnewer ${COMPIL_DIR}/$1 -name $2 -print ) 
     149      [ ${#mycpp} -ne 0 ] && cpeval ${myarch} ${COMPIL_DIR}/$1 
    141150       fi 
    142        # has the cpp keys file been updated since we copied the arch file in ${COMPIL_DIR}? 
    143        mycpp=$( find -L ${COMPIL_DIR} -cnewer ${COMPIL_DIR}/$1 -name $2 -print ) 
    144        [ ${#mycpp} -ne 0 ] && cpeval ${myarch} ${COMPIL_DIR}/$1 
    145151       # has myarch file been updated since we copied it in ${COMPIL_DIR}? 
    146152       myarch=$( find -L ${MAIN_DIR}/ARCH -cnewer ${COMPIL_DIR}/$1 -name arch-${3}.fcm -print ) 
    … …  
    150156   fi 
    151157    else 
    152    ls -l ${COMPIL_DIR}/$2 | sed -e "s/.* -> //" > ${COMPIL_DIR}/cpp.history 
     158   if [ "$2" != "nocpp" ]  
     159   then 
     160       ls -l ${COMPIL_DIR}/$2 | sed -e "s/.* -> //" > ${COMPIL_DIR}/cpp.history 
     161   fi 
    153162   echo ${myarch} > ${COMPIL_DIR}/arch.history 
    154163   cpeval ${myarch} ${COMPIL_DIR}/$1 
    … …  
    157166 
    158167#- do we need xios library? 
    159 use_iom=$( sed -e "s/#.*$//" ${COMPIL_DIR}/$2 | grep -c key_iomput ) 
     168if [ "$2" != "nocpp" ]  
     169then 
     170    use_iom=$( sed -e "s/#.*$//" ${COMPIL_DIR}/$2 | grep -c key_iomput ) 
     171else 
     172    use_iom=0 
     173fi 
    160174have_lxios=$( sed -e "s/#.*$//" ${COMPIL_DIR}/$1 | grep -c "\-lxios" ) 
    161175if [[ ( $use_iom -eq 0 ) && ( $have_lxios -ge 1 ) ]] 
    … …  
    166180 
    167181#- do we need oasis libraries? 
    168 use_oasis=$( sed -e "s/#.*$//" ${COMPIL_DIR}/$2 | grep -c key_oasis3 ) 
     182if [ "$2" != "nocpp" ]  
     183then 
     184    use_oasis=$( sed -e "s/#.*$//" ${COMPIL_DIR}/$2 | grep -c key_oasis3 ) 
     185else 
     186    use_oasis=0 
     187fi 
    169188for liboa in psmile.MPI1 mct mpeu scrip mpp_io 
    170189do 

  • branches/2013/dev_UKMO_2013/NEMOGCM/TOOLS/MISCELLANEOUS/chk_iomput.sh

    r3974 r4171  
    5959#------------------------------------------------ 
    6060# 
    61 external=$( grep -c "<field_definition.* src=" $xmlfile ) 
     61external=$( grep -c "<field_definition  *\([^ ].* \)*src=" $xmlfile ) 
    6262if [ $external -eq 1 ] 
    6363then 
    64     xmlfield_def=$( grep "<field_definition.* src=" $xmlfile | sed -e 's/.*src="\([^"]*\)".*/\1/' ) 
     64    xmlfield_def=$( grep "<field_definition  *\([^ ].* \)*src=" $xmlfile | sed -e 's/.*src="\([^"]*\)".*/\1/' ) 
    6565    xmlfield_def=$( dirname $xmlfile )/$xmlfield_def    
    6666else 
    6767    xmlfield_def=$xmlfile 
    6868fi 
    69 [ $inxml -eq 1 ] && grep "< *field * id *=" $xmlfield_def 
     69[ $inxml -eq 1 ] && grep "< *field  *\([^ ].* \)*id *=" $xmlfield_def 
    7070[ $insrc -eq 1 ] && find $srcdir -name "*.[Ffh]90" -exec grep -iH "^[^\!]*call  *iom_put *(" {} \; 
    7171[ $(( $insrc + $inxml )) -ge 1 ] && exit 
    … …  
    9191# list of variables defined in the xml file 
    9292# 
    93 varlistxml=$( grep "< *field.* id *=" $xmlfield_def  | sed -e "s/^.*< *field.* id *= *[\"\']\([^\"\']*\)[\"\'].*/\1/" | sort -d ) 
     93varlistxml=$( grep "< *field  *\([^ ].* \)*id *=" $xmlfield_def  | sed -e "s/^.*< *field .*id *= *[\"\']\([^\"\']*\)[\"\'].*/\1/" | sort -d ) 
    9494# 
    9595# list of variables to be outputed in the xml file 
    9696# 
    97 varlistout=$( grep "< *field.* field_ref *=" $xmlfile  | sed -e "s/^.*< *field.* field_ref *= *[\"\']\([^\"\']*\)[\"\'].*/\1/" | sort -d ) 
     97varlistout=$( grep "< *field  *\([^ ].* \)*field_ref *=" $xmlfile  | sed -e "s/^.*< *field .*field_ref *= *[\"\']\([^\"\']*\)[\"\'].*/\1/" | sort -d ) 
    9898# 
    9999echo "--------------------------------------------------" 

  • branches/2013/dev_UKMO_2013/NEMOGCM/TOOLS/maketools

    r3294 r4171  
    146146 
    147147#- When used for the first time, choose a compiler --- 
    148 . ${COMPIL_DIR}/Fcheck_archfile.sh arch_tools.fcm ${CMP_NAM} || exit 
     148. ${COMPIL_DIR}/Fcheck_archfile.sh arch_tools.fcm nocpp ${CMP_NAM} || exit 
    149149 
    150150#- Choose a default tool if needed --- 
Note: See TracChangeset for help on using the changeset viewer.