Changeset 3370 for branches/2012/dev_r3337_NOCS10_ICB/NEMOGCM/NEMO/OPA_SRC
- Timestamp:
- 2012-04-30T10:27:44+02:00 (12 years ago)
- Location:
- branches/2012/dev_r3337_NOCS10_ICB/NEMOGCM/NEMO/OPA_SRC
- Files:
-
- 12 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2012/dev_r3337_NOCS10_ICB/NEMOGCM/NEMO/OPA_SRC/ICB/icb_oce.F90
r3359 r3370 3 3 !!====================================================================== 4 4 !! *** MODULE icb_oce *** 5 !! Ocean physics: declare variables for iceberg tracking5 !! Icebergs: declare variables for iceberg tracking 6 6 !!====================================================================== 7 !! History : 3.3.1 ! 2010-01 (Martin&Adcroft) Original code 8 !! - ! 2011-03 (Madec) Part conversion to NEMO form 9 !! - ! Removal of mapping from another grid 10 !! - ! 2011-04 (Alderson) Extensive rewrite 11 !! - ! Split into separate modules 7 !! History : 3.3 ! 2010-01 (T. Martin & A. Adcroft) Original code 8 !! - ! 2011-03 (G. Madec) Part conversion to NEMO form 9 !! - ! Removal of mapping from another grid 10 !! - ! 2011-04 (S. Alderson) Extensive rewrite ; Split into separate modules 12 11 !!---------------------------------------------------------------------- 13 12 !! … … 16 15 !! 17 16 !! A single iceberg is held as an instance of type 'iceberg' 18 !! This type defines a linked list, so each instance contains a pointer to the previous and next19 !! icebergs in the list17 !! This type defines a linked list, so each instance contains a pointer 18 !! to the previous and next icebergs in the list 20 19 !! 21 20 !! Type 'icebergs' is a convenience container for all relevant arrays 22 21 !! It contains one pointer to an 'iceberg' instance representing all icebergs in the processor 23 22 !! 24 !! Each iceberg has a position represented as a real cartesian coordinate which is fractional25 !! grid cell, centred on T points; so an iceberg position of (1.0,1.0) lies exactly on the first26 !! T point and the Tcell spans 0.5 to 1.5 in each direction23 !! Each iceberg has a position represented as a real cartesian coordinate which is 24 !! fractional grid cell, centred on T-points; so an iceberg position of (1.0,1.0) lies 25 !! exactly on the first T-point and the T-cell spans 0.5 to 1.5 in each direction 27 26 !! 28 27 !! Each iceberg is assigned a unique id even in MPI 29 !! This consists of an array of integers: the first element is used to label, the second and30 !! subsequent elements are used to count the number of times the first element wraps around31 !! a ll possible values within the valid size for this datatype.28 !! This consists of an array of integers: the first element is used to label, the second 29 !! and subsequent elements are used to count the number of times the first element wraps 30 !! around all possible values within the valid size for this datatype. 32 31 !! Labelling is done by starting the first label in each processor (even when only one) 33 32 !! as narea, and then incrementing by jpnij (i.e. the total number of processors. … … 58 57 END TYPE icebergs_gridded 59 58 60 !! point type is just a convenience storage type for properties of an individual iceberg 61 TYPE, PUBLIC :: point 62 INTEGER :: year 63 REAL(wp) :: xi, yj ! iceberg coordinates in the (i,j) referential (global) 64 REAL(wp) :: e1, e2 ! horizontal scale factors at the iceberg position 65 REAL(wp) :: lon, lat, day 66 REAL(wp) :: mass, thickness, width, length, uvel, vvel 67 REAL(wp) :: uo, vo, ui, vi, ua, va, ssh_x, ssh_y, sst, cn, hi 68 REAL(wp) :: mass_of_bits, heat_density 59 TYPE, PUBLIC :: point !: properties of an individual iceberg (position, mass, size, etc...) 60 INTEGER :: year 61 REAL(wp) :: xi , yj ! iceberg coordinates in the (i,j) referential (global) 62 REAL(wp) :: e1 , e2 ! horizontal scale factors at the iceberg position 63 REAL(wp) :: lon, lat, day ! geographic position 64 REAL(wp) :: mass, thickness, width, length, uvel, vvel ! iceberg physical properties 65 REAL(wp) :: uo, vo, ui, vi, ua, va, ssh_x, ssh_y, sst, cn, hi ! properties of iceberg environment 66 REAL(wp) :: mass_of_bits, heat_density 69 67 END TYPE point 70 68 71 !! iceberg type is really just a linked list 72 TYPE, PUBLIC :: iceberg 73 ! pointers to previous and next unique icebergs in linked list 74 TYPE(iceberg), POINTER :: prev=>NULL(), next=>NULL() 75 ! variables which do not change for this iceberg 76 INTEGER, DIMENSION(nkounts) :: number 77 REAL(wp) :: mass_scaling 78 ! variables which change with time are held in a separate type 79 TYPE(point), POINTER :: current_point => NULL() 69 TYPE, PUBLIC :: iceberg !: linked list defining all the icebergs present in the model domain 70 TYPE(iceberg), POINTER :: prev=>NULL(), next=>NULL() ! pointers to previous and next unique icebergs in linked list 71 INTEGER, DIMENSION(nkounts) :: number ! variables which do not change for this iceberg 72 REAL(wp) :: mass_scaling ! - - - - 73 TYPE(point), POINTER :: current_point => NULL() ! variables which change with time are held in a separate type 80 74 END TYPE iceberg 81 75 82 !! master instances of gridded and linked list iceberg types83 TYPE(icebergs_gridded), POINTER :: berg_grid84 TYPE(iceberg) , POINTER :: first_berg => NULL()85 76 86 !! parameters controlling iceberg characteristics and modelling 87 REAL(wp) :: berg_dt ! Time-step between iceberg CALLs (should make adaptive?) 88 INTEGER :: current_year 89 REAL(wp) :: current_yearday ! 1.00-365.99 90 REAL(wp), DIMENSION(:), POINTER :: first_width, first_length 91 LOGICAL :: l_restarted_bergs=.FALSE. ! Indicate whether we read state from a restart or not 92 ! ! arbitrary numbers for diawri entry 77 TYPE(icebergs_gridded), POINTER :: berg_grid !: master instance of gridded iceberg type 78 TYPE(iceberg) , POINTER :: first_berg => NULL() !: master instance of linked list iceberg type 79 80 ! !!! parameters controlling iceberg characteristics and modelling 81 REAL(wp) :: berg_dt !: Time-step between iceberg CALLs (should make adaptive?) 82 INTEGER :: current_year !: 83 REAL(wp) :: current_yearday !: 1.00-365.99 84 REAL(wp), DIMENSION(:), POINTER :: first_width, first_length !: 85 LOGICAL :: l_restarted_bergs=.FALSE. ! Indicate whether we read state from a restart or not 86 ! ! arbitrary numbers for diawri entry 93 87 REAL(wp), DIMENSION(nclasses), PUBLIC :: class_num=(/ 1, 2, 3, 4, 5, 6, 7, 8, 9, 10 /) 94 88 95 89 ! Extra arrays with bigger halo, needed when interpolating forcing onto iceberg position 96 90 ! particularly for MPP when iceberg can lie inside T grid but outside U, V, or f grid 97 REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE :: uo_e 98 REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE :: vo_e 99 REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE :: ff_e 100 REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE :: ua_e 101 REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE :: va_e 91 REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE :: uo_e, vo_e 92 REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE :: ff_e 93 REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE :: ua_e, va_e 102 94 #if defined key_lim2 || defined key_lim3 103 REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE :: ui_e 104 REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE :: vi_e 105 REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE :: ssh_e 95 REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE :: ui_e, vi_e 96 REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE :: ssh_e 106 97 #endif 107 98 108 99 !!gm almost all those PARAM ARE defined in NEMO 109 REAL(wp), PUBLIC, PARAMETER :: pp_rho_ice = 916.7_wp ! Density of fresh ice @ 0oC (kg/m^3)110 REAL(wp), PUBLIC, PARAMETER :: pp_rho_water = 999.8_wp ! Density of fresh water @ 0oC (kg/m^3)111 REAL(wp), PUBLIC, PARAMETER :: pp_rho_air = 1.1_wp ! Density of air @ 0oC (kg/m^3) ???112 REAL(wp), PUBLIC, PARAMETER :: pp_rho_seawater = 1025._wp ! Approx. density of surface sea water @ 0oC (kg/m^3)100 REAL(wp), PUBLIC, PARAMETER :: pp_rho_ice = 916.7_wp !: Density of fresh ice @ 0oC [kg/m^3] 101 REAL(wp), PUBLIC, PARAMETER :: pp_rho_water = 999.8_wp !: Density of fresh water @ 0oC [kg/m^3] 102 REAL(wp), PUBLIC, PARAMETER :: pp_rho_air = 1.1_wp !: Density of air @ 0oC [kg/m^3] 103 REAL(wp), PUBLIC, PARAMETER :: pp_rho_seawater = 1025._wp !: Approx. density of surface sea water @ 0oC [kg/m^3] 113 104 !!gm end 114 REAL(wp), PUBLIC, PARAMETER :: pp_Cd_av = 1.3_wp ! (Vertical) Drag coefficient between bergs and atmos (?)115 REAL(wp), PUBLIC, PARAMETER :: pp_Cd_ah = 0.0055_wp ! (Horizontal) Drag coefficient between bergs and atmos (?)116 REAL(wp), PUBLIC, PARAMETER :: pp_Cd_wv = 0.9_wp ! (Vertical) Drag coefficient between bergs and ocean (?)117 REAL(wp), PUBLIC, PARAMETER :: pp_Cd_wh = 0.0012_wp ! (Horizontal) Drag coefficient between bergs and ocean (?)118 REAL(wp), PUBLIC, PARAMETER :: pp_Cd_iv = 0.9_wp ! (Vertical) Drag coefficient between bergs and sea-ice (?)119 !TOM> no horizontal drag for sea ice! real, PARAMETER :: pp_Cd_ih=0.0012 ! (Horizontal) Drag coefficient between bergs and sea-ice (?) 105 REAL(wp), PUBLIC, PARAMETER :: pp_Cd_av = 1.3_wp !: (Vertical) Drag coefficient between bergs and atmos 106 REAL(wp), PUBLIC, PARAMETER :: pp_Cd_ah = 0.0055_wp !: (lateral ) Drag coefficient between bergs and atmos 107 REAL(wp), PUBLIC, PARAMETER :: pp_Cd_wv = 0.9_wp !: (Vertical) Drag coefficient between bergs and ocean 108 REAL(wp), PUBLIC, PARAMETER :: pp_Cd_wh = 0.0012_wp !: (lateral ) Drag coefficient between bergs and ocean 109 REAL(wp), PUBLIC, PARAMETER :: pp_Cd_iv = 0.9_wp !: (Vertical) Drag coefficient between bergs and sea-ice 110 !TOM> no horizontal drag for sea ice! real, PARAMETER :: pp_Cd_ih=0.0012 ! (lateral) Drag coeff. between bergs and sea-ice 120 111 121 ! !** namberg namelist parameters (and defaults) ** 122 LOGICAL, PUBLIC :: ln_bergdia = .true. ! Calculate budgets 123 INTEGER, PUBLIC :: nn_verbose_level = 0 ! Turn on debugging when level > 0 124 INTEGER, PUBLIC :: nn_test_icebergs = 0 ! Create icebergs in absence of a restart file from the supplied class number 125 REAL(wp), PUBLIC, DIMENSION(4) :: rn_test_box = (/ 0._wp, 1._wp, 0._wp, 1._wp /) 126 ! lon1,lon2,lat1,lat2 box to create them in 127 INTEGER, PUBLIC :: nn_sample_rate = 0 ! Timesteps between sampling of position for trajectory storage 128 INTEGER, PUBLIC :: nn_verbose_write = 15 ! timesteps between verbose messages 129 REAL(wp), PUBLIC :: rn_rho_bergs = 850._wp ! Density of icebergs 130 REAL(wp), PUBLIC :: rn_LoW_ratio = 1.5_wp ! Initial ratio L/W for newly calved icebergs 131 REAL(wp), PUBLIC :: rn_bits_erosion_fraction = 0. ! Fraction of erosion melt flux to divert to bergy bits 132 REAL(wp), PUBLIC :: rn_sicn_shift = 0._wp ! Shift of sea-ice concentration in erosion flux modulation (0<sicn_shift<1) 133 LOGICAL, PUBLIC :: ln_operator_splitting = .true. ! Use first order operator splitting for thermodynamics 134 LOGICAL, PUBLIC :: ln_passive_mode = .false. ! iceberg - ocean decoupling 135 LOGICAL, PUBLIC :: ln_time_average_weight = .false. ! Time average the weight on the ocean !!gm I don't understand that ! 136 REAL(wp), PUBLIC :: rn_speed_limit = 0._wp ! CFL speed limit for a berg 112 ! !!* namberg namelist parameters (and defaults) ** 113 LOGICAL , PUBLIC :: ln_bergdia = .true. !: Calculate budgets 114 INTEGER , PUBLIC :: nn_verbose_level = 0 !: Turn on debugging when level > 0 115 INTEGER , PUBLIC :: nn_test_icebergs = 0 !: Create icebergs in absence of a restart file from the supplied class nb 116 REAL(wp), PUBLIC, DIMENSION(4) :: rn_test_box = (/ 0._wp, 1._wp, 0._wp, 1._wp /) !: lon1,lon2,lat1,lat2 box to create them in 117 INTEGER , PUBLIC :: nn_sample_rate = 0 !: Timesteps between sampling of position for trajectory storage 118 INTEGER , PUBLIC :: nn_verbose_write = 15 !: timesteps between verbose messages 119 REAL(wp), PUBLIC :: rn_rho_bergs = 850._wp !: Density of icebergs 120 REAL(wp), PUBLIC :: rn_LoW_ratio = 1.5_wp !: Initial ratio L/W for newly calved icebergs 121 REAL(wp), PUBLIC :: rn_bits_erosion_fraction = 0. !: Fraction of erosion melt flux to divert to bergy bits 122 REAL(wp), PUBLIC :: rn_sicn_shift = 0._wp !: Shift of sea-ice concentration in erosion flux modulation (0<sicn_shift<1) 123 LOGICAL , PUBLIC :: ln_operator_splitting = .true. !: Use first order operator splitting for thermodynamics 124 LOGICAL , PUBLIC :: ln_passive_mode = .false. !: iceberg - ocean decoupling 125 LOGICAL , PUBLIC :: ln_time_average_weight = .false. !: Time average the weight on the ocean !!gm I don't understand that ! 126 REAL(wp), PUBLIC :: rn_speed_limit = 0._wp !: CFL speed limit for a berg 137 127 ! 138 ! ! Mass thresholds between iceberg classes (kg) 139 REAL(wp), DIMENSION(nclasses), PUBLIC :: rn_initial_mass=(/8.8e7, 4.1e8, 3.3e9, 1.8e10, 3.8e10, 7.5e10, 1.2e11, & 140 2.2e11, 3.9e11, 7.4e11/) 141 ! ! Fraction of calving to apply to this class (non-dim) 142 REAL(wp), DIMENSION(nclasses), PUBLIC :: rn_distribution=(/0.25, 0.12, 0.15, 0.18, 0.12, 0.07, 0.03, 0.03, 0.03, 0.02/) 143 ! ! Ratio between effective and real iceberg mass (non-dim) 144 REAL(wp), DIMENSION(nclasses), PUBLIC :: rn_mass_scaling=(/2000, 200, 50, 20, 10, 5, 2, 1, 1, 1/) 145 ! ! Total thickness of newly calved bergs (m) 146 REAL(wp), DIMENSION(nclasses), PUBLIC :: rn_initial_thickness=(/40., 67., 133., 175., 250., 250., 250., 250., 250., 250./) 128 ! ! Mass thresholds between iceberg classes [kg] 129 REAL(wp), DIMENSION(nclasses), PUBLIC :: rn_initial_mass = (/8.8e07, 4.1e08, 3.3e09, 1.8e10, 3.8e10, & 130 & 7.5e10, 1.2e11, 2.2e11, 3.9e11, 7.4e11/) 131 ! ! Fraction of calving to apply to this class [non-dim] 132 REAL(wp), DIMENSION(nclasses), PUBLIC :: rn_distribution = (/0.25 , 0.12 , 0.15 , 0.18 , 0.12 , & 133 & 0.07 , 0.03 , 0.03 , 0.03 , 0.02 /) 134 ! ! Ratio between effective and real iceberg mass (non-dim) 135 REAL(wp), DIMENSION(nclasses), PUBLIC :: rn_mass_scaling = (/2000. , 200. , 50. , 20. , 10. , & 136 & 5. , 2. , 1. , 1. , 1. /) 137 ! ! Total thickness of newly calved bergs [m] 138 REAL(wp), DIMENSION(nclasses), PUBLIC :: rn_initial_thickness = (/ 40. , 67. , 133. , 175. , 250. , & 139 & 250. , 250. , 250. , 250. , 250./) 147 140 148 141 ! Single instance of an icebergs type initialised in icebergs_init and updated in icebergs_run 149 REAL(wp), ALLOCATABLE, PUBLIC, SAVE, DIMENSION(:,:) :: src_calving, src_calving_hflx150 INTEGER , PUBLIC, SAVE :: numicb151 INTEGER , PUBLIC, SAVE, DIMENSION(nkounts) :: num_bergs152 INTEGER , PUBLIC, SAVE :: nicbdi, nicbei, nicbdj, nicbej153 REAL(wp) , PUBLIC, SAVE :: ricb_left, ricb_right154 INTEGER , PUBLIC, SAVE :: nicbpack155 INTEGER , PUBLIC, SAVE :: nktberg, nknberg156 INTEGER , ALLOCATABLE , PUBLIC, SAVE, DIMENSION(:) :: nicbfldpts157 INTEGER , ALLOCATABLE , PUBLIC, SAVE, DIMENSION(:) :: nicbflddest158 INTEGER , ALLOCATABLE , PUBLIC, SAVE, DIMENSION(:) :: nicbfldproc142 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: src_calving, src_calving_hflx !: ??? 143 INTEGER , PUBLIC , SAVE :: numicb !: ??? 144 INTEGER , PUBLIC , SAVE, DIMENSION(nkounts) :: num_bergs !: ??? 145 INTEGER , PUBLIC , SAVE :: nicbdi, nicbei, nicbdj, nicbej !: ??? 146 REAL(wp), PUBLIC , SAVE :: ricb_left, ricb_right !: ??? 147 INTEGER , PUBLIC , SAVE :: nicbpack !: ??? 148 INTEGER , PUBLIC , SAVE :: nktberg, nknberg !: ??? 149 INTEGER , PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: nicbfldpts !: ??? 150 INTEGER , PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: nicbflddest !: ??? 151 INTEGER , PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: nicbfldproc !: ??? 159 152 160 TYPE(FLD), ALLOCATABLE, PUBLIC , DIMENSION(:) :: sf_icb !structure: file information, fields read153 TYPE(FLD), ALLOCATABLE, PUBLIC , DIMENSION(:) :: sf_icb !: structure: file information, fields read 161 154 162 155 END MODULE icb_oce -
branches/2012/dev_r3337_NOCS10_ICB/NEMOGCM/NEMO/OPA_SRC/ICB/icbclv.F90
r3359 r3370 3 3 !!====================================================================== 4 4 !! *** MODULE icbclv *** 5 !! Ocean physics: calving routines for iceberg calving5 !! Icebergs: calving routines for iceberg calving 6 6 !!====================================================================== 7 7 !! History : 3.3.1 ! 2010-01 (Martin&Adcroft) Original code … … 12 12 !!---------------------------------------------------------------------- 13 13 !!---------------------------------------------------------------------- 14 !! accumulate_calving 15 !! icb_gen : generate test icebergs14 !! accumulate_calving : 15 !! icb_gen : generate test icebergs 16 16 !!---------------------------------------------------------------------- 17 17 USE par_oce ! NEMO parameters … … 21 21 USE lbclnk ! NEMO boundary exchanges for gridded data 22 22 23 USE icb_oce ! define iceberg arrays24 USE icbdia ! iceberg utility routines23 USE icb_oce ! iceberg variables 24 USE icbdia ! iceberg diagnostics 25 25 USE icbutl ! iceberg utility routines 26 26 … … 31 31 PUBLIC calve_icebergs ! routine called in icbrun.F90 module 32 32 33 !!---------------------------------------------------------------------- 34 !! NEMO/OPA 3.3 , NEMO Consortium (2011) 35 !! $Id:$ 36 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 37 !!---------------------------------------------------------------------- 33 38 CONTAINS 34 39 … … 41 46 !! ** input : - ? 42 47 !!---------------------------------------------------------------------- 43 !44 48 INTEGER, INTENT(in) :: kt 45 49 ! … … 63 67 !do jn=1, nclasses 64 68 ! where (berg_grid%calving==0.) berg_grid%stored_ice(:,:,jn)=0. 65 !end do69 !end do 66 70 DO jj = 2, jpjm1 67 71 DO ji = 2, jpim1 … … 109 113 !! is considered 110 114 !!---------------------------------------------------------------------- 115 INTEGER :: ji, jj, jn ! dummy loop indices 116 INTEGER :: icnt, icntmax 117 TYPE(iceberg) :: newberg 118 TYPE(point) :: newpt 119 REAL(wp) :: zddt, zcalved_to_berg, zheat_to_berg 120 !!---------------------------------------------------------------------- 111 121 ! 112 INTEGER :: icnt,icntmax 113 TYPE(iceberg) :: newberg 114 TYPE(point) :: newpt 115 INTEGER :: ji, jj, jn ! dummy loop indices 116 REAL(wp) :: zddt, zcalved_to_berg, zheat_to_berg 117 !!---------------------------------------------------------------------- 118 119 icntmax = 0 120 122 icntmax = 0 123 ! 121 124 DO jn = 1, nclasses 122 125 DO jj = nicbdj, nicbej … … 128 131 DO WHILE (berg_grid%stored_ice(ji,jj,jn) >= rn_initial_mass(jn) * rn_mass_scaling(jn) ) 129 132 ! 130 newpt%lon = glamt(ji,jj) ! at t-point (centre of the cell)133 newpt%lon = glamt(ji,jj) ! at t-point (centre of the cell) 131 134 newpt%lat = gphit(ji,jj) 132 135 newpt%xi = REAL( nimpp+ji-1, wp ) 133 136 newpt%yj = REAL( njmpp+jj-1, wp ) 134 137 ! 135 newpt%uvel = 0._wp ! initially at rest138 newpt%uvel = 0._wp ! initially at rest 136 139 newpt%vvel = 0._wp 137 140 ! ! set berg characteristics 138 newpt%mass = rn_initial_mass (jn)139 newpt%thickness = rn_initial_thickness(jn)140 newpt%width = first_width(jn)141 newpt%length = first_length(jn)141 newpt%mass = rn_initial_mass (jn) 142 newpt%thickness = rn_initial_thickness(jn) 143 newpt%width = first_width (jn) 144 newpt%length = first_length (jn) 142 145 newberg%mass_scaling = rn_mass_scaling (jn) 143 newpt%mass_of_bits = 0._wp ! no bergy146 newpt%mass_of_bits = 0._wp ! no bergy 144 147 ! 145 148 newpt%year = current_year … … 157 160 berg_grid%stored_heat(ji,jj) = berg_grid%stored_heat(ji,jj) - zheat_to_berg 158 161 ! ! Stored mass 159 berg_grid%stored_ice(ji,jj,jn) 162 berg_grid%stored_ice(ji,jj,jn) = berg_grid%stored_ice(ji,jj,jn) - zcalved_to_berg 160 163 ! 161 164 zddt = zddt + berg_dt * 2._wp / 17._wp ! Minor offset to start day !!gm why??? … … 171 174 DO jn = 1,nclasses 172 175 CALL lbc_lnk( berg_grid%stored_ice(:,:,jn), 'T', 1._wp ) 173 END DO176 END DO 174 177 CALL lbc_lnk( berg_grid%stored_heat, 'T', 1._wp ) 175 178 ! … … 178 181 END SUBROUTINE calve_icebergs 179 182 183 !!====================================================================== 180 184 END MODULE icbclv -
branches/2012/dev_r3337_NOCS10_ICB/NEMOGCM/NEMO/OPA_SRC/ICB/icbdia.F90
r3361 r3370 3 3 !!====================================================================== 4 4 !! *** MODULE icbdia *** 5 !! Ocean physics: initialise variables for iceberg budgets and diagnostics5 !! Icebergs: initialise variables for iceberg budgets and diagnostics 6 6 !!====================================================================== 7 !! History : 3.3 .1 ! 2010-01 (Martin&Adcroft) Original code8 !! - 9 !! - 10 !! - 11 !! - 12 !! - 13 !! - 7 !! History : 3.3 ! 2010-01 (Martin, Adcroft) Original code 8 !! - ! 2011-03 (Madec) Part conversion to NEMO form 9 !! - ! Removal of mapping from another grid 10 !! - ! 2011-04 (Alderson) Split into separate modules 11 !! - ! 2011-05 (Alderson) Budgets are now all here with lots 12 !! - ! of silly routines to call to get values in 13 !! - ! from the right points in the code 14 14 !!---------------------------------------------------------------------- 15 15 !!---------------------------------------------------------------------- … … 17 17 !! icb_budget_init : initialise iceberg budgeting 18 18 !!---------------------------------------------------------------------- 19 USE par_oce ! nemoparameters19 USE par_oce ! ocean parameters 20 20 USE dom_oce ! ocean domain 21 21 USE in_out_manager ! nemo IO 22 USE lib_mpp 23 USE iom 24 25 USE icb_oce ! define iceberg arrays 22 USE lib_mpp ! MPP library 23 USE iom ! I/O library 24 USE icb_oce ! iceberg variables 26 25 USE icbutl ! iceberg utility routines 27 26 … … 29 28 PRIVATE 30 29 31 REAL(wp), DIMENSION(:,:) , POINTER, PUBLIC :: berg_melt=>NULL() ! Melting+erosion rate of icebergs (kg/s/m^2) 32 REAL(wp), DIMENSION(:,:) , POINTER, PUBLIC :: buoy_melt=>NULL() ! Buoyancy component of melting rate (kg/s/m^2) 33 REAL(wp), DIMENSION(:,:) , POINTER, PUBLIC :: eros_melt=>NULL() ! Erosion component of melting rate (kg/s/m^2) 34 REAL(wp), DIMENSION(:,:) , POINTER, PUBLIC :: conv_melt=>NULL() ! Convective component of melting rate (kg/s/m^2) 35 REAL(wp), DIMENSION(:,:) , POINTER, PUBLIC :: bits_src=>NULL() ! Mass flux from berg erosion into bergy bits (kg/s/m^2) 36 REAL(wp), DIMENSION(:,:) , POINTER, PUBLIC :: bits_melt=>NULL() ! Melting rate of bergy bits (kg/s/m^2) 37 REAL(wp), DIMENSION(:,:) , POINTER, PUBLIC :: bits_mass=>NULL() ! Mass distribution of bergy bits (kg/s/m^2) 38 REAL(wp), DIMENSION(:,:) , POINTER, PUBLIC :: virtual_area=>NULL() ! Virtual surface coverage by icebergs (m^2) 39 REAL(wp), DIMENSION(:,:) , POINTER, PUBLIC :: berg_mass=>NULL() ! Mass distribution (kg/m^2) 40 REAL(wp), DIMENSION(:,:,:), POINTER, PUBLIC :: real_calving=>NULL() ! Calving rate into iceberg class at calving locations (kg/s) 41 REAL(wp), DIMENSION(:,:) , POINTER, PRIVATE :: tmpc=>NULL() ! Temporary work space 42 REAL(wp), DIMENSION(:) , POINTER, PRIVATE :: rsumbuf=>NULL() ! Temporary work space to reduce mpp exchanges 43 INTEGER , DIMENSION(:) , POINTER, PRIVATE :: nsumbuf=>NULL() ! Temporary work space to reduce mpp exchanges 44 45 REAL(wp) , PRIVATE :: berg_melt_net 46 REAL(wp) , PRIVATE :: bits_src_net 47 REAL(wp) , PRIVATE :: bits_melt_net 48 REAL(wp) , PRIVATE :: bits_mass_start, bits_mass_end 49 REAL(wp) , PRIVATE :: floating_heat_start, floating_heat_end 50 REAL(wp) , PRIVATE :: floating_mass_start, floating_mass_end 51 REAL(wp) , PRIVATE :: bergs_mass_start, bergs_mass_end 52 REAL(wp) , PRIVATE :: stored_start, stored_heat_start 53 REAL(wp) , PRIVATE :: stored_end , stored_heat_end 54 REAL(wp) , PRIVATE :: calving_src_net, calving_out_net 55 REAL(wp) , PRIVATE :: calving_src_heat_net, calving_out_heat_net 56 REAL(wp) , PRIVATE :: calving_src_heat_used_net 57 REAL(wp) , PRIVATE :: calving_rcv_net, calving_ret_net, calving_used_net 58 REAL(wp) , PRIVATE :: heat_to_bergs_net, heat_to_ocean_net, melt_net 59 REAL(wp) , PRIVATE :: calving_to_bergs_net 60 61 INTEGER , PRIVATE :: nbergs_start, nbergs_end, nbergs_calved 62 INTEGER , PRIVATE :: nbergs_melted 63 INTEGER , PRIVATE :: nspeeding_tickets 64 INTEGER , DIMENSION(nclasses) , PRIVATE :: nbergs_calved_by_class 65 66 PUBLIC icb_budget_end ! routine called in icbrun.F90 module 67 PUBLIC icb_budget_init ! routine called in icbini.F90 module 68 PUBLIC icb_budget ! routine called in icbrun.F90 module 69 PUBLIC icb_budget_step ! routine called in icbrun.F90 module 70 PUBLIC icb_budget_put ! routine called in icbrun.F90 module 71 PUBLIC melt_budget ! routine called in icbthm.F90 module 72 PUBLIC size_budget ! routine called in icbthm.F90 module 73 PUBLIC speed_budget ! routine called in icbdyn.F90 module 74 PUBLIC calving_budget ! routine called in icbclv.F90 module 75 PUBLIC incoming_budget ! routine called in icbclv.F90 module 76 30 PUBLIC icb_budget_end ! routine called in icbrun.F90 module 31 PUBLIC icb_budget_init ! routine called in icbini.F90 module 32 PUBLIC icb_budget ! routine called in icbrun.F90 module 33 PUBLIC icb_budget_step ! routine called in icbrun.F90 module 34 PUBLIC icb_budget_put ! routine called in icbrun.F90 module 35 PUBLIC melt_budget ! routine called in icbthm.F90 module 36 PUBLIC size_budget ! routine called in icbthm.F90 module 37 PUBLIC speed_budget ! routine called in icbdyn.F90 module 38 PUBLIC calving_budget ! routine called in icbclv.F90 module 39 PUBLIC incoming_budget ! routine called in icbclv.F90 module 40 41 REAL(wp), DIMENSION(:,:) , POINTER, PUBLIC :: berg_melt => NULL() ! Melting+erosion rate of icebergs [kg/s/m2] 42 REAL(wp), DIMENSION(:,:) , POINTER, PUBLIC :: buoy_melt => NULL() ! Buoyancy component of melting rate [kg/s/m2] 43 REAL(wp), DIMENSION(:,:) , POINTER, PUBLIC :: eros_melt => NULL() ! Erosion component of melting rate [kg/s/m2] 44 REAL(wp), DIMENSION(:,:) , POINTER, PUBLIC :: conv_melt => NULL() ! Convective component of melting rate [kg/s/m2] 45 REAL(wp), DIMENSION(:,:) , POINTER, PUBLIC :: bits_src => NULL() ! Mass flux from berg erosion into bergy bits [kg/s/m2] 46 REAL(wp), DIMENSION(:,:) , POINTER, PUBLIC :: bits_melt => NULL() ! Melting rate of bergy bits [kg/s/m2] 47 REAL(wp), DIMENSION(:,:) , POINTER, PUBLIC :: bits_mass => NULL() ! Mass distribution of bergy bits [kg/s/m2] 48 REAL(wp), DIMENSION(:,:) , POINTER, PUBLIC :: virtual_area => NULL() ! Virtual surface coverage by icebergs [m2] 49 REAL(wp), DIMENSION(:,:) , POINTER, PUBLIC :: berg_mass => NULL() ! Mass distribution [kg/m2] 50 REAL(wp), DIMENSION(:,:,:), POINTER, PUBLIC :: real_calving => NULL() ! Calving rate into iceberg class at 51 ! ! calving locations [kg/s] 52 53 REAL(wp), DIMENSION(:,:) , POINTER :: tmpc => NULL() ! Temporary work space 54 REAL(wp), DIMENSION(:) , POINTER :: rsumbuf => NULL() ! Temporary work space to reduce mpp exchanges 55 INTEGER , DIMENSION(:) , POINTER :: nsumbuf => NULL() ! Temporary work space to reduce mpp exchanges 56 57 REAL(wp) :: berg_melt_net 58 REAL(wp) :: bits_src_net 59 REAL(wp) :: bits_melt_net 60 REAL(wp) :: bits_mass_start , bits_mass_end 61 REAL(wp) :: floating_heat_start , floating_heat_end 62 REAL(wp) :: floating_mass_start , floating_mass_end 63 REAL(wp) :: bergs_mass_start , bergs_mass_end 64 REAL(wp) :: stored_start , stored_heat_start 65 REAL(wp) :: stored_end , stored_heat_end 66 REAL(wp) :: calving_src_net , calving_out_net 67 REAL(wp) :: calving_src_heat_net, calving_out_heat_net 68 REAL(wp) :: calving_src_heat_used_net 69 REAL(wp) :: calving_rcv_net , calving_ret_net , calving_used_net 70 REAL(wp) :: heat_to_bergs_net, heat_to_ocean_net, melt_net 71 REAL(wp) :: calving_to_bergs_net 72 73 INTEGER :: nbergs_start, nbergs_end, nbergs_calved 74 INTEGER :: nbergs_melted 75 INTEGER :: nspeeding_tickets 76 INTEGER , DIMENSION(nclasses) :: nbergs_calved_by_class 77 78 !!---------------------------------------------------------------------- 79 !! NEMO/OPA 3.3 , NEMO Consortium (2011) 80 !! $Id:$ 81 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 82 !!---------------------------------------------------------------------- 77 83 CONTAINS 78 84 79 !!-------------------------------------------------------------------------80 81 85 SUBROUTINE icb_budget_end 82 83 IF( .NOT. ln_bergdia ) RETURN 86 !!---------------------------------------------------------------------- 87 ! 88 IF( .NOT. ln_bergdia ) RETURN 84 89 DEALLOCATE( berg_melt ) 85 90 DEALLOCATE( buoy_melt ) … … 97 102 DEALLOCATE( nsumbuf ) 98 103 ENDIF 99 104 ! 100 105 END SUBROUTINE icb_budget_end 101 106 … … 103 108 104 109 SUBROUTINE icb_budget_init( ) 105 110 !!---------------------------------------------------------------------- 111 !!---------------------------------------------------------------------- 112 ! 106 113 IF( .NOT. ln_bergdia ) RETURN 107 ALLOCATE( berg_melt (jpi,jpj) ) ; berg_melt (:,:)= 0._wp108 ALLOCATE( buoy_melt (jpi,jpj) ) ; buoy_melt (:,:)= 0._wp109 ALLOCATE( eros_melt (jpi,jpj) ) ; eros_melt (:,:)= 0._wp110 ALLOCATE( conv_melt (jpi,jpj) ) ; conv_melt (:,:)= 0._wp111 ALLOCATE( bits_src (jpi,jpj) ) ; bits_src (:,:) = 0._wp112 ALLOCATE( bits_melt (jpi,jpj) ) ; bits_melt (:,:) = 0._wp113 ALLOCATE( bits_mass (jpi,jpj) ) ; bits_mass (:,:) = 0._wp114 ALLOCATE( virtual_area (jpi,jpj) ) ; virtual_area (:,:)= 0._wp115 ALLOCATE( berg_mass (jpi,jpj) ) ; berg_mass (:,:)= 0._wp116 ALLOCATE( real_calving (jpi,jpj,nclasses) ) ; real_calving (:,:,:)=0.117 ALLOCATE( tmpc(jpi,jpj) ) ; tmpc (:,:)=0.114 ALLOCATE( berg_melt (jpi,jpj) ) ; berg_melt (:,:) = 0._wp 115 ALLOCATE( buoy_melt (jpi,jpj) ) ; buoy_melt (:,:) = 0._wp 116 ALLOCATE( eros_melt (jpi,jpj) ) ; eros_melt (:,:) = 0._wp 117 ALLOCATE( conv_melt (jpi,jpj) ) ; conv_melt (:,:) = 0._wp 118 ALLOCATE( bits_src (jpi,jpj) ) ; bits_src (:,:) = 0._wp 119 ALLOCATE( bits_melt (jpi,jpj) ) ; bits_melt (:,:) = 0._wp 120 ALLOCATE( bits_mass (jpi,jpj) ) ; bits_mass (:,:) = 0._wp 121 ALLOCATE( virtual_area (jpi,jpj) ) ; virtual_area(:,:) = 0._wp 122 ALLOCATE( berg_mass (jpi,jpj) ) ; berg_mass (:,:) = 0._wp 123 ALLOCATE( real_calving (jpi,jpj,nclasses) ) ; real_calving(:,:,:) = 0._wp 124 ALLOCATE( tmpc(jpi,jpj) ) ; tmpc (:,:) = 0._wp 118 125 119 126 nbergs_start = 0 … … 167 174 bits_mass_start = rsumbuf(3) 168 175 ENDIF 169 176 ! 170 177 END SUBROUTINE icb_budget_init 171 178 172 !!-------------------------------------------------------------------------173 179 174 180 SUBROUTINE icb_budget( ld_budge ) 175 ! Arguments 176 LOGICAL, INTENT(in) :: ld_budge 177 ! Local variables 178 INTEGER :: ik 179 REAL(wp) :: zunused_calving, ztmpsum, zgrdd_berg_mass, zgrdd_bits_mass 180 181 IF( .NOT. ln_bergdia ) RETURN 182 183 zunused_calving = SUM( berg_grid%calving(:,:) ) 184 ztmpsum = SUM( berg_grid%floating_melt(:,:) * e1e2t(:,:) * tmask_i(:,:) ) 185 melt_net = melt_net + ztmpsum * berg_dt 186 calving_out_net = calving_out_net + ( zunused_calving + ztmpsum ) * berg_dt 187 ztmpsum = SUM( berg_melt(:,:) * e1e2t(:,:) * tmask_i(:,:) ) 188 berg_melt_net = berg_melt_net + ztmpsum * berg_dt 189 ztmpsum = SUM( bits_src(:,:) * e1e2t(:,:) * tmask_i(:,:) ) 190 bits_src_net = bits_src_net + ztmpsum * berg_dt 191 ztmpsum = SUM( bits_melt(:,:) * e1e2t(:,:) * tmask_i(:,:) ) 192 bits_melt_net = bits_melt_net + ztmpsum * berg_dt 193 ztmpsum = SUM( src_calving(:,:) * tmask_i(:,:) ) 194 calving_ret_net = calving_ret_net + ztmpsum * berg_dt 195 ztmpsum = SUM( berg_grid%calving_hflx(:,:) * e1e2t(:,:) * tmask_i(:,:) ) 196 calving_out_heat_net = calving_out_heat_net + ztmpsum * berg_dt ! Units of J 181 !!---------------------------------------------------------------------- 182 !!---------------------------------------------------------------------- 183 LOGICAL, INTENT(in) :: ld_budge 184 ! 185 INTEGER :: ik 186 REAL(wp) :: zunused_calving, ztmpsum, zgrdd_berg_mass, zgrdd_bits_mass 187 !!---------------------------------------------------------------------- 188 ! 189 IF( .NOT. ln_bergdia ) RETURN 190 191 zunused_calving = SUM( berg_grid%calving(:,:) ) 192 ztmpsum = SUM( berg_grid%floating_melt(:,:) * e1e2t(:,:) * tmask_i(:,:) ) 193 melt_net = melt_net + ztmpsum * berg_dt 194 calving_out_net = calving_out_net + ( zunused_calving + ztmpsum ) * berg_dt 195 ztmpsum = SUM( berg_melt(:,:) * e1e2t(:,:) * tmask_i(:,:) ) 196 berg_melt_net = berg_melt_net + ztmpsum * berg_dt 197 ztmpsum = SUM( bits_src(:,:) * e1e2t(:,:) * tmask_i(:,:) ) 198 bits_src_net = bits_src_net + ztmpsum * berg_dt 199 ztmpsum = SUM( bits_melt(:,:) * e1e2t(:,:) * tmask_i(:,:) ) 200 bits_melt_net = bits_melt_net + ztmpsum * berg_dt 201 ztmpsum = SUM( src_calving(:,:) * tmask_i(:,:) ) 202 calving_ret_net = calving_ret_net + ztmpsum * berg_dt 203 ztmpsum = SUM( berg_grid%calving_hflx(:,:) * e1e2t(:,:) * tmask_i(:,:) ) 204 calving_out_heat_net = calving_out_heat_net + ztmpsum * berg_dt ! Units of J 197 205 198 206 IF( ld_budge ) THEN 199 stored_end 200 stored_heat_end 201 floating_mass_end 202 bergs_mass_end 203 bits_mass_end 204 floating_heat_end 205 206 nbergs_end 207 zgrdd_berg_mass 208 zgrdd_bits_mass 207 stored_end = SUM( berg_grid%stored_ice(:,:,:) ) 208 stored_heat_end = SUM( berg_grid%stored_heat(:,:) ) 209 floating_mass_end = sum_mass( first_berg ) 210 bergs_mass_end = sum_mass( first_berg,justbergs=.true. ) 211 bits_mass_end = sum_mass( first_berg,justbits=.true. ) 212 floating_heat_end = sum_heat( first_berg ) 213 214 nbergs_end = count_bergs() 215 zgrdd_berg_mass = SUM( berg_mass (:,:)*e1e2t(:,:)*tmask_i(:,:) ) 216 zgrdd_bits_mass = SUM( bits_mass(:,:)*e1e2t(:,:)*tmask_i(:,:) ) 209 217 210 218 IF( lk_mpp ) THEN … … 263 271 nsumbuf(3) = nbergs_melted 264 272 nsumbuf(4) = nspeeding_tickets 265 DO ik = 1, nclasses273 DO ik = 1, nclasses 266 274 nsumbuf(4+ik) = nbergs_calved_by_class(ik) 267 END DO275 END DO 268 276 269 277 CALL mpp_sum( nsumbuf(1:nclasses+4), nclasses+4 ) … … 355 363 bits_src_net = 0._wp 356 364 ENDIF 357 365 ! 358 366 END SUBROUTINE icb_budget 359 367 360 !!-------------------------------------------------------------------------361 368 362 369 SUBROUTINE icb_budget_step 370 !!---------------------------------------------------------------------- 363 371 !! things to reset at the beginning of each timestep 364 372 !! this probably screws up fields going to diawri, so needs to be looked at - sga 365 373 !!---------------------------------------------------------------------- 374 ! 366 375 IF( .NOT. ln_bergdia ) RETURN 367 376 berg_melt (:,:) = 0._wp … … 375 384 virtual_area (:,:) = 0._wp 376 385 real_calving (:,:,:) = 0._wp 377 386 ! 378 387 END SUBROUTINE icb_budget_step 379 388 380 !!-------------------------------------------------------------------------381 389 382 390 SUBROUTINE icb_budget_put 383 384 IF( .NOT. ln_bergdia ) RETURN 385 CALL iom_put( "berg_melt" , berg_melt (:,:) ) ! 'Melt rate of icebergs' , 'kg/m2/s' 386 CALL iom_put( "berg_buoy_melt" , buoy_melt (:,:) ) ! 'Buoyancy component of iceberg melt rate' , 'kg/m2/s' 387 CALL iom_put( "berg_eros_melt" , eros_melt (:,:) ) ! 'Erosion component of iceberg melt rate' , 'kg/m2/s' 388 CALL iom_put( "berg_conv_melt" , conv_melt (:,:) ) ! 'Convective component of iceberg melt rate', 'kg/m2/s' 389 CALL iom_put( "berg_virtual_area" , virtual_area (:,:) ) ! 'Virtual coverage by icebergs' , 'm2' 390 CALL iom_put( "bits_src" , bits_src (:,:) ) ! 'Mass source of bergy bits' , 'kg/m2/s' 391 CALL iom_put( "bits_melt" , bits_melt (:,:) ) ! 'Melt rate of bergy bits' , 'kg/m2/s' 392 CALL iom_put( "bits_mass" , bits_mass (:,:) ) ! 'Bergy bit density field' , 'kg/m2' 393 CALL iom_put( "berg_mass" , berg_mass (:,:) ) ! 'Iceberg density field' , 'kg/m2' 394 CALL iom_put( "berg_real_calving" , real_calving (:,:,:) ) ! 'Calving into iceberg class' , 'kg/s' 395 391 !!---------------------------------------------------------------------- 392 !!---------------------------------------------------------------------- 393 ! 394 IF( .NOT. ln_bergdia ) RETURN !!gm useless iom will control whether it is output or not 395 ! 396 CALL iom_put( "berg_melt" , berg_melt (:,:) ) ! Melt rate of icebergs [kg/m2/s] 397 CALL iom_put( "berg_buoy_melt" , buoy_melt (:,:) ) ! Buoyancy component of iceberg melt rate [kg/m2/s] 398 CALL iom_put( "berg_eros_melt" , eros_melt (:,:) ) ! Erosion component of iceberg melt rate [kg/m2/s] 399 CALL iom_put( "berg_conv_melt" , conv_melt (:,:) ) ! Convective component of iceberg melt rate [kg/m2/s] 400 CALL iom_put( "berg_virtual_area", virtual_area(:,:) ) ! Virtual coverage by icebergs [m2] 401 CALL iom_put( "bits_src" , bits_src (:,:) ) ! Mass source of bergy bits [kg/m2/s] 402 CALL iom_put( "bits_melt" , bits_melt (:,:) ) ! Melt rate of bergy bits [kg/m2/s] 403 CALL iom_put( "bits_mass" , bits_mass (:,:) ) ! Bergy bit density field [kg/m2] 404 CALL iom_put( "berg_mass" , berg_mass (:,:) ) ! Iceberg density field [kg/m2] 405 CALL iom_put( "berg_real_calving", real_calving(:,:,:) ) ! Calving into iceberg class [kg/s] 406 ! 396 407 END SUBROUTINE icb_budget_put 397 408 398 !!-------------------------------------------------------------------------399 409 400 410 SUBROUTINE calving_budget( ki, kj, kn, pcalved, pheated ) 411 !!---------------------------------------------------------------------- 412 !!---------------------------------------------------------------------- 401 413 INTEGER, INTENT(in) :: ki, kj, kn 402 414 REAL(wp), INTENT(in) :: pcalved 403 415 REAL(wp), INTENT(in) :: pheated 404 416 !!---------------------------------------------------------------------- 417 ! 405 418 IF( .NOT. ln_bergdia ) RETURN 406 419 real_calving(ki,kj,kn) = real_calving(ki,kj,kn) + pcalved / berg_dt … … 409 422 calving_to_bergs_net = calving_to_bergs_net + pcalved 410 423 heat_to_bergs_net = heat_to_bergs_net + pheated 411 424 ! 412 425 END SUBROUTINE calving_budget 413 426 414 !!-------------------------------------------------------------------------415 427 416 428 SUBROUTINE incoming_budget( kt, pcalving_used, pheat_used ) 429 !!---------------------------------------------------------------------- 430 !!---------------------------------------------------------------------- 417 431 INTEGER , INTENT(in) :: kt 418 432 REAL(wp), INTENT(in) :: pcalving_used 419 433 REAL(wp), DIMENSION(:,:), INTENT(in) :: pheat_used 420 434 !!---------------------------------------------------------------------- 435 ! 421 436 IF( .NOT. ln_bergdia ) RETURN 422 437 ! 423 438 IF( kt == nit000 ) THEN 424 439 stored_start = SUM( berg_grid%stored_ice(:,:,:) ) 425 440 IF( lk_mpp ) CALL mpp_sum( stored_start ) 426 441 WRITE(numicb,'(a,es13.6,a)') 'accumulate_calving: initial stored mass=',stored_start,' kg' 427 442 ! 428 443 stored_heat_start = SUM( berg_grid%stored_heat(:,:) ) 429 444 IF( lk_mpp ) CALL mpp_sum( stored_heat_start ) 430 445 WRITE(numicb,'(a,es13.6,a)') 'accumulate_calving: initial stored heat=',stored_heat_start,' J' 431 446 ENDIF 432 447 ! 433 448 calving_rcv_net = calving_rcv_net + SUM( berg_grid%calving(:,:) ) * berg_dt 434 449 calving_src_net = calving_rcv_net 435 450 calving_src_heat_net = calving_src_heat_net + & 436 451 & SUM( berg_grid%calving_hflx(:,:) * e1e2t(:,:) ) * berg_dt ! Units of J 437 452 calving_used_net = calving_used_net + pcalving_used * berg_dt 438 453 calving_src_heat_used_net = calving_src_heat_used_net + SUM( pheat_used(:,:) ) 439 454 ! 440 455 END SUBROUTINE incoming_budget 441 456 442 !!-------------------------------------------------------------------------443 457 444 458 SUBROUTINE size_budget(ki, kj, pWn, pLn, pAbits, & 445 459 & pmass_scale, pMnew, pnMbits, pz1_e1e2) 446 INTEGER, INTENT(in) :: ki, kj 447 REAL(wp), INTENT(in) :: pWn, pLn, pAbits, pmass_scale, pMnew, pnMbits, pz1_e1e2 448 460 !!---------------------------------------------------------------------- 461 !!---------------------------------------------------------------------- 462 INTEGER, INTENT(in) :: ki, kj 463 REAL(wp), INTENT(in) :: pWn, pLn, pAbits, pmass_scale, pMnew, pnMbits, pz1_e1e2 464 !!---------------------------------------------------------------------- 465 ! 449 466 IF( .NOT. ln_bergdia ) RETURN 450 467 virtual_area(ki,kj) = virtual_area(ki,kj) + ( pWn * pLn + pAbits ) * pmass_scale ! m^2 451 468 berg_mass(ki,kj) = berg_mass(ki,kj) + pMnew * pz1_e1e2 ! kg/m2 452 469 bits_mass(ki,kj) = bits_mass(ki,kj) + pnMbits * pz1_e1e2 ! kg/m2 453 470 ! 454 471 END SUBROUTINE size_budget 455 472 456 !!-------------------------------------------------------------------------457 473 458 474 SUBROUTINE speed_budget() 459 475 !!---------------------------------------------------------------------- 476 !!---------------------------------------------------------------------- 477 ! 460 478 IF( .NOT. ln_bergdia ) RETURN 461 479 nspeeding_tickets = nspeeding_tickets + 1 462 480 ! 463 481 END SUBROUTINE speed_budget 464 482 465 !!-------------------------------------------------------------------------466 483 467 484 SUBROUTINE melt_budget(ki, kj, pmnew, pheat, pmass_scale, & 468 485 & pdM, pdMbitsE, pdMbitsM, pdMb, pdMe, & 469 & pdMv, pz1_dt_e1e2) 470 471 INTEGER :: ki, kj 472 REAL(wp), INTENT(in) :: pmnew, pheat, pmass_scale 473 REAL(wp), INTENT(in) :: pdM, pdMbitsE, pdMbitsM, pdMb, pdMe, pdMv, pz1_dt_e1e2 474 486 & pdMv, pz1_dt_e1e2 ) 487 !!---------------------------------------------------------------------- 488 !!---------------------------------------------------------------------- 489 INTEGER , INTENT(in) :: ki, kj 490 REAL(wp), INTENT(in) :: pmnew, pheat, pmass_scale 491 REAL(wp), INTENT(in) :: pdM, pdMbitsE, pdMbitsM, pdMb, pdMe, pdMv, pz1_dt_e1e2 492 !!---------------------------------------------------------------------- 493 ! 475 494 IF( .NOT. ln_bergdia ) RETURN 476 495 … … 483 502 heat_to_ocean_net = heat_to_ocean_net + pheat * pmass_scale * berg_dt ! J 484 503 IF( pmnew <= 0._wp ) nbergs_melted = nbergs_melted + 1 ! Delete the berg if completely melted 485 504 ! 486 505 END SUBROUTINE melt_budget 487 506 488 !!-------------------------------------------------------------------------489 507 490 508 SUBROUTINE report_state( cd_budgetstr, cd_budgetunits, cd_startstr, pstartval, cd_endstr, & 491 509 & pendval, cd_delstr, kbergs ) 492 ! Arguments 510 !!---------------------------------------------------------------------- 511 !!---------------------------------------------------------------------- 493 512 CHARACTER*(*), INTENT(in) :: cd_budgetstr, cd_budgetunits, cd_startstr, cd_endstr, cd_delstr 494 513 REAL(wp), INTENT(in) :: pstartval, pendval 495 514 INTEGER, INTENT(in), OPTIONAL :: kbergs 496 515 !!---------------------------------------------------------------------- 516 ! 497 517 IF ( PRESENT(kbergs) ) THEN 498 518 WRITE(numicb,100) cd_budgetstr // ' state:', & … … 510 530 END SUBROUTINE report_state 511 531 512 !!-------------------------------------------------------------------------513 532 514 533 SUBROUTINE report_consistant( cd_budgetstr, cd_budgetunits, cd_startstr, pstartval, cd_endstr, pendval) 515 ! Arguments 534 !!---------------------------------------------------------------------- 535 !!---------------------------------------------------------------------- 516 536 CHARACTER*(*), INTENT(in) :: cd_budgetstr, cd_budgetunits, cd_startstr, cd_endstr 517 537 REAL(wp), INTENT(in) :: pstartval, pendval 518 538 !!---------------------------------------------------------------------- 539 ! 519 540 WRITE(numicb,200) cd_budgetstr // ' check:', & 520 541 cd_startstr, pstartval, cd_budgetunits, & … … 524 545 END SUBROUTINE report_consistant 525 546 526 !!-------------------------------------------------------------------------527 547 528 548 SUBROUTINE report_budget( cd_budgetstr, cd_budgetunits, cd_instr, pinval, cd_outstr, & 529 549 & poutval, cd_delstr, pstartval, pendval) 530 ! Arguments 550 !!---------------------------------------------------------------------- 551 !!---------------------------------------------------------------------- 531 552 CHARACTER*(*), INTENT(in) :: cd_budgetstr, cd_budgetunits, cd_instr, cd_outstr, cd_delstr 532 553 REAL(wp), INTENT(in) :: pinval, poutval, pstartval, pendval 533 ! Local variables554 ! 534 555 REAL(wp) :: zval 535 556 !!---------------------------------------------------------------------- 557 ! 536 558 zval = ( ( pendval - pstartval ) - ( pinval - poutval ) ) / & 537 559 & MAX( 1.e-30, MAX( abs( pendval - pstartval ) , ABS( pinval - poutval ) ) ) 538 560 539 561 WRITE(numicb,200) cd_budgetstr // ' budget:', & 540 cd_instr // ' in', pinval, cd_budgetunits, & 541 cd_outstr // ' out', poutval, cd_budgetunits, & 542 'Delta ' // cd_delstr, pinval-poutval, cd_budgetunits, & 543 'error', zval, 'nd' 544 200 FORMAT(a19,3(a18,"=",es14.7,x,a2,:,","),a8,"=",es10.3,x,a2) 562 & cd_instr // ' in', pinval, cd_budgetunits, & 563 & cd_outstr // ' out', poutval, cd_budgetunits, & 564 & 'Delta ' // cd_delstr, pinval-poutval, cd_budgetunits, & 565 & 'error', zval, 'nd' 566 200 FORMAT(a19,3(a18,"=",es14.7,x,a2,:,","),a8,"=",es10.3,x,a2) 567 ! 545 568 END SUBROUTINE report_budget 546 569 547 !!-------------------------------------------------------------------------548 570 549 571 SUBROUTINE report_istate( cd_budgetstr, cd_startstr, pstartval, cd_endstr, pendval, cd_delstr) 550 ! Arguments 572 !!---------------------------------------------------------------------- 573 !!---------------------------------------------------------------------- 551 574 CHARACTER*(*), INTENT(in) :: cd_budgetstr, cd_startstr, cd_endstr, cd_delstr 552 575 INTEGER, INTENT(in) :: pstartval, pendval 553 576 ! 554 577 WRITE(numicb,100) cd_budgetstr // ' state:', & 555 cd_startstr // ' start', pstartval, & 556 cd_endstr // ' end', pendval, & 557 cd_delstr // 'Delta', pendval-pstartval 558 100 FORMAT(a19,3(a18,"=",i14,x,:,",")) 578 & cd_startstr // ' start', pstartval, & 579 & cd_endstr // ' end', pendval, & 580 & cd_delstr // 'Delta', pendval-pstartval 581 100 FORMAT(a19,3(a18,"=",i14,x,:,",")) 582 ! 559 583 END SUBROUTINE report_istate 560 584 561 !!-------------------------------------------------------------------------562 585 563 586 SUBROUTINE report_ibudget( cd_budgetstr, cd_instr, pinval, cd_outstr, poutval, & 564 587 & cd_delstr, pstartval, pendval) 565 ! Arguments 588 !!---------------------------------------------------------------------- 589 !!---------------------------------------------------------------------- 566 590 CHARACTER*(*), INTENT(in) :: cd_budgetstr, cd_instr, cd_outstr, cd_delstr 567 591 INTEGER, INTENT(in) :: pinval, poutval, pstartval, pendval 568 592 !!---------------------------------------------------------------------- 593 ! 569 594 WRITE(numicb,200) cd_budgetstr // ' budget:', & 570 595 cd_instr // ' in', pinval, & … … 574 599 200 FORMAT(a19,10(a18,"=",i14,x,:,",")) 575 600 END SUBROUTINE report_ibudget 576 !!------------------------------------------------------------------------- 577 601 602 !!====================================================================== 578 603 END MODULE icbdia -
branches/2012/dev_r3337_NOCS10_ICB/NEMOGCM/NEMO/OPA_SRC/ICB/icbdyn.F90
r3359 r3370 3 3 !!====================================================================== 4 4 !! *** MODULE icbdyn *** 5 !! Ocean physics: time stepping routine for iceberg tracking5 !! Iceberg: time stepping routine for iceberg tracking 6 6 !!====================================================================== 7 7 !! History : 3.3.1 ! 2010-01 (Martin&Adcroft) Original code … … 16 16 USE dom_oce ! NEMO ocean domain 17 17 USE phycst ! NEMO physical constants 18 18 ! 19 19 USE icb_oce ! define iceberg arrays 20 20 USE icbutl ! iceberg utility routines … … 26 26 PUBLIC evolve_icebergs ! routine called in icbrun.F90 module 27 27 28 !!---------------------------------------------------------------------- 29 !! NEMO/OPA 3.3 , NEMO Consortium (2011) 30 !! $Id:$ 31 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 32 !!---------------------------------------------------------------------- 28 33 CONTAINS 29 34 … … 36 41 !! ** Method : - blah blah 37 42 !!---------------------------------------------------------------------- 38 !39 43 REAL(wp) :: zuvel1 , zvvel1 , zu1, zv1, zax1, zay1, zxi1 , zyj1 40 44 REAL(wp) :: zuvel2 , zvvel2 , zu2, zv2, zax2, zay2, zxi2 , zyj2 41 45 REAL(wp) :: zuvel3 , zvvel3 , zu3, zv3, zax3, zay3, zxi3 , zyj3 42 46 REAL(wp) :: zuvel4 , zvvel4 , zu4, zv4, zax4, zay4, zxi4 , zyj4 43 REAL(wp) :: zuvel_n, zvvel_n , zxi_n, zyj_n47 REAL(wp) :: zuvel_n, zvvel_n, zxi_n , zyj_n 44 48 REAL(wp) :: zdt, zdt_2, zdt_6, ze1, ze2 45 49 LOGICAL :: ll_bounced … … 48 52 !!---------------------------------------------------------------------- 49 53 50 ! 4th order Runge-Kutta to solve: 51 ! d/dt X = V, d/dt V = A 52 ! with I.C.'s: 53 ! X=X1 and V=V1 54 ! 4th order Runge-Kutta to solve: d/dt X = V, d/dt V = A 55 ! with I.C.'s: X=X1 and V=V1 54 56 ! 55 57 ! ; A1=A(X1,V1) … … 82 84 83 85 ! !** A1 = A(X1,V1) 84 CALL accel( 85 & berg ,zyj1, ze2, zvvel1, zvvel1, zay1, zdt_2 )86 CALL accel( berg , zxi1, ze1, zuvel1, zuvel1, zax1, & 87 & zyj1, ze2, zvvel1, zvvel1, zay1, zdt_2 ) 86 88 ! 87 89 zu1 = zuvel1 / ze1 !** V1 in d(i,j)/dt … … 99 101 100 102 ! !** A2 = A(X2,V2) 101 CALL accel( 102 & berg ,zyj2, ze2, zvvel2, zvvel1, zay2, zdt_2 )103 CALL accel( berg , zxi2, ze1, zuvel2, zuvel1, zax2, & 104 & zyj2, ze2, zvvel2, zvvel1, zay2, zdt_2 ) 103 105 ! 104 106 zu2 = zuvel2 / ze1 !** V2 in d(i,j)/dt … … 115 117 116 118 ! !** A3 = A(X3,V3) 117 CALL accel( 118 & berg ,zyj3, ze2, zvvel3, zvvel1, zay3, zdt )119 CALL accel( berg , zxi3, ze1, zuvel3, zuvel1, zax3, & 120 & zyj3, ze2, zvvel3, zvvel1, zay3, zdt ) 119 121 ! 120 122 zu3 = zuvel3 / ze1 !** V3 in d(i,j)/dt … … 131 133 132 134 ! !** A4 = A(X4,V4) 133 CALL accel( 134 & berg ,zyj4, ze2, zvvel4, zvvel1, zay4, zdt )135 CALL accel( berg , zxi4, ze1, zuvel4, zuvel1, zax4, & 136 & zyj4, ze2, zvvel4, zvvel1, zay4, zdt ) 135 137 136 138 zu4 = zuvel4 / ze1 !** V4 in d(i,j)/dt … … 174 176 !! ** Method : - blah blah 175 177 !!---------------------------------------------------------------------- 176 REAL(wp), INTENT(inout) :: pi , pj! current iceberg position177 REAL(wp), INTENT(in ) :: pi0, pj0! previous iceberg position178 REAL(wp), INTENT(inout) :: pu , pv! current iceberg velocities179 LOGICAL , INTENT( out) :: ld_bounced! bounced indicator178 REAL(wp), INTENT(inout) :: pi , pj ! current iceberg position 179 REAL(wp), INTENT(in ) :: pi0, pj0 ! previous iceberg position 180 REAL(wp), INTENT(inout) :: pu , pv ! current iceberg velocities 181 LOGICAL , INTENT( out) :: ld_bounced ! bounced indicator 180 182 ! 181 183 INTEGER :: ii, ii0 … … 231 233 232 234 233 SUBROUTINE accel( 234 & berg ,pyj, pe2, pvvel, pvvel0, pay, pdt )235 SUBROUTINE accel( berg , pxi, pe1, puvel, puvel0, pax, & 236 & pyj, pe2, pvvel, pvvel0, pay, pdt ) 235 237 !!---------------------------------------------------------------------- 236 238 !! *** ROUTINE accel *** … … 240 242 !! ** Method : - blah blah 241 243 !!---------------------------------------------------------------------- 242 TYPE(iceberg ), POINTER ::berg243 REAL(wp) , INTENT(in ):: pxi , pyj ! berg position in (i,j) referential244 REAL(wp) , INTENT(in ):: puvel , pvvel ! berg velocity [m/s]245 REAL(wp) , INTENT(in ):: puvel0, pvvel0 ! initial berg velocity [m/s]246 REAL(wp) , INTENT( out):: pe1, pe2 ! horizontal scale factor at (xi,yj)247 REAL(wp) , INTENT(inout):: pax, pay ! berg acceleration248 REAL(wp) , INTENT(in ) :: pdt! berg time step244 TYPE(iceberg ), POINTER, INTENT(in ) :: berg ! berg 245 REAL(wp) , INTENT(in ) :: pxi , pyj ! berg position in (i,j) referential 246 REAL(wp) , INTENT(in ) :: puvel , pvvel ! berg velocity [m/s] 247 REAL(wp) , INTENT(in ) :: puvel0, pvvel0 ! initial berg velocity [m/s] 248 REAL(wp) , INTENT( out) :: pe1, pe2 ! horizontal scale factor at (xi,yj) 249 REAL(wp) , INTENT(inout) :: pax, pay ! berg acceleration 250 REAL(wp) , INTENT(in ) :: pdt ! berg time step 249 251 ! 250 252 REAL(wp), PARAMETER :: pp_alpha = 0._wp ! … … 270 272 271 273 zM = berg%current_point%mass 272 zT = berg%current_point%thickness 273 zD = ( rn_rho_bergs / pp_rho_seawater ) * zT ! draught (keel depth)274 zT = berg%current_point%thickness ! total thickness 275 zD = ( rn_rho_bergs / pp_rho_seawater ) * zT ! draught (keel depth) 274 276 zF = zT - zD ! freeboard 275 277 zW = berg%current_point%width … … 281 283 ! Wave radiation 282 284 zuwave = zua - zuo ; zvwave = zva - zvo ! Use wind speed rel. to ocean for wave model 283 zwmod = zuwave*zuwave + zvwave*zvwave ! The wave amplitude and length depend on the current;284 ! ! wind speed relative to the ocean. Actually wmod is wmod**2 here.285 zampl = 0.5 * 0.02025 * zwmod ! This is "a", the wave amplitude286 zLwavelength = 0.32 * zwmod ! Surface wave length fitted to data in table at287 ! ! http://www4.ncsu.edu/eos/users/c/ceknowle/public/chapter10/part2.html285 zwmod = zuwave*zuwave + zvwave*zvwave ! The wave amplitude and length depend on the current; 286 ! ! wind speed relative to the ocean. Actually wmod is wmod**2 here. 287 zampl = 0.5 * 0.02025 * zwmod ! This is "a", the wave amplitude 288 zLwavelength = 0.32 * zwmod ! Surface wave length fitted to data in table at 289 ! ! http://www4.ncsu.edu/eos/users/c/ceknowle/public/chapter10/part2.html 288 290 zLcutoff = 0.125 * zLwavelength 289 291 zLtop = 0.25 * zLwavelength 290 292 zCr = pp_Cr0 * MIN( MAX( 0., (zL-zLcutoff) / ((zLtop-zLcutoff)+1.e-30)) , 1.) ! Wave radiation coefficient 291 ! ! fitted to graph from Carrieres et al., POAC Drift Model.293 ! ! fitted to graph from Carrieres et al., POAC Drift Model. 292 294 zwave_rad = 0.5 * pp_rho_seawater / zM * zCr * grav * zampl * MIN( zampl,zF ) * (2.*zW*zL) / (zW+zL) 293 zwmod = SQRT( zua*zua + zva*zva ) ! Wind speed295 zwmod = SQRT( zua*zua + zva*zva ) ! Wind speed 294 296 IF( zwmod /= 0._wp ) THEN 295 zuwave = zua/zwmod ! Wave radiation force acts in wind direction ... !!gm this should be the wind rel. to ocean ?297 zuwave = zua/zwmod ! Wave radiation force acts in wind direction ... !!gm this should be the wind rel. to ocean ? 296 298 zvwave = zva/zwmod 297 299 ELSE … … 325 327 zaxe = zaxe + zff*pvvel0 326 328 zaye = zaye - zff*puvel0 327 else328 zaxe =zaxe+zff*pvvel329 zaye =zaye-zff*puvel330 endif331 if (pp_beta>0.) then! If implicit, use time-level (n) rather than RK4 latest332 zaxe =zaxe-zdrag_ocn*(puvel0-zuo) -zdrag_atm*(puvel0-zua) -zdrag_ice*(puvel0-zui)333 zaye =zaye-zdrag_ocn*(pvvel0-zvo) -zdrag_atm*(pvvel0-zva) -zdrag_ice*(pvvel0-zvi)334 else335 zaxe =zaxe-zdrag_ocn*(puvel-zuo) -zdrag_atm*(puvel-zua) -zdrag_ice*(puvel-zui)336 zaye =zaye-zdrag_ocn*(pvvel-zvo) -zdrag_atm*(pvvel-zva) -zdrag_ice*(pvvel-zvi)329 ELSE 330 zaxe = zaxe + zff*pvvel 331 zaye = zaye - zff*puvel 332 ENDIF 333 IF( pp_beta > 0._wp ) THEN ! If implicit, use time-level (n) rather than RK4 latest 334 zaxe = zaxe - zdrag_ocn*(puvel0-zuo) - zdrag_atm*(puvel0-zua) -zdrag_ice*(puvel0-zui) 335 zaye = zaye - zdrag_ocn*(pvvel0-zvo) - zdrag_atm*(pvvel0-zva) -zdrag_ice*(pvvel0-zvi) 336 ELSE 337 zaxe = zaxe - zdrag_ocn*(puvel -zuo) - zdrag_atm*(puvel -zua) -zdrag_ice*(puvel -zui) 338 zaye = zaye - zdrag_ocn*(pvvel -zvo) - zdrag_atm*(pvvel -zva) -zdrag_ice*(pvvel -zvi) 337 339 endif 338 340 … … 340 342 IF( pp_alpha + pp_beta > 0._wp ) THEN 341 343 zlambda = zdrag_ocn + zdrag_atm + zdrag_ice 342 zA11 = 1.+pp_beta*pdt*zlambda343 zA12 =pp_alpha*pdt*zff344 zdetA = 1._wp / ( zA11*zA11 + zA12*zA12 )345 pax = zdetA * ( zA11*zaxe+zA12*zaye )346 pay = zdetA * ( zA11*zaye-zA12*zaxe )344 zA11 = 1._wp + pp_beta *pdt*zlambda 345 zA12 = pp_alpha*pdt*zff 346 zdetA = 1._wp / ( zA11*zA11 + zA12*zA12 ) 347 pax = zdetA * ( zA11*zaxe + zA12*zaye ) 348 pay = zdetA * ( zA11*zaye - zA12*zaxe ) 347 349 else 348 pax =zaxe ; pay=zaye350 pax = zaxe ; pay = zaye 349 351 endif 350 352 … … 354 356 END DO ! itloop 355 357 356 IF( rn_speed_limit > 0. ) THEN ! Limit speed of bergs based on a CFL criteria (if asked)358 IF( rn_speed_limit > 0._wp ) THEN ! Limit speed of bergs based on a CFL criteria (if asked) 357 359 zspeed = SQRT( zuveln*zuveln + zvveln*zvveln ) ! Speed of berg 358 IF( zspeed > 0. ) THEN360 IF( zspeed > 0._wp ) THEN 359 361 zloc_dx = MIN( pe1, pe2 ) ! minimum grid spacing 360 362 zspeed_new = zloc_dx / pdt * rn_speed_limit ! Speed limit as a factor of dx / dt … … 374 376 END SUBROUTINE accel 375 377 378 !!====================================================================== 376 379 END MODULE icbdyn -
branches/2012/dev_r3337_NOCS10_ICB/NEMOGCM/NEMO/OPA_SRC/ICB/icbini.F90
r3369 r3370 3 3 !!====================================================================== 4 4 !! *** MODULE icbini *** 5 !! Ocean physics: initialise variables for iceberg tracking5 !! Icebergs: initialise variables for iceberg tracking 6 6 !!====================================================================== 7 !! History : 3.3.1 ! 2010-01 (Martin&Adcroft) Original code 8 !! - ! 2011-03 (Madec) Part conversion to NEMO form 9 !! - ! Removal of mapping from another grid 10 !! - ! 2011-04 (Alderson) Split into separate modules 11 !! - ! Restore restart routines 12 !! - ! 2011-05 (Alderson) generate_test_icebergs restored 13 !! - ! 2011-05 (Alderson) new forcing arrays with extra halo 14 !! - ! 2011-05 (Alderson) north fold exchange arrays added 7 !! History : - ! 2010-01 (T. Martin & A. Adcroft) Original code 8 !! 3.3 ! 2011-03 (G. Madec) Part conversion to NEMO form ; Removal of mapping from another grid 9 !! - ! 2011-04 (S. Alderson) Split into separate modules ; Restore restart routines 10 !! - ! 2011-05 (S. Alderson) generate_test_icebergs restored ; new forcing arrays with extra halo ; 11 !! - ! north fold exchange arrays added 15 12 !!---------------------------------------------------------------------- 16 13 !!---------------------------------------------------------------------- 17 !! icb_init : initialise 14 !! icb_init : initialise icebergs 18 15 !! icb_gen : generate test icebergs 19 16 !! icb_nam : read iceberg namelist … … 22 19 USE in_out_manager ! IO routines and numout in particular 23 20 USE lib_mpp ! mpi library and lk_mpp in particular 24 USE sbc_oce 25 USE iom 26 USE fldread 27 USE lbclnk 28 21 USE sbc_oce ! ocean: surface boundary condition 22 USE iom ! IOM library 23 USE fldread ! field read 24 USE lbclnk ! lateral boundary condition - MPP link 25 ! 29 26 USE icb_oce ! define iceberg arrays 30 27 USE icbutl ! iceberg utility routines … … 36 33 PRIVATE 37 34 38 CHARACTER(len=100), PRIVATE :: cn_dir = './' ! Root directory for location of icb files 39 TYPE(FLD_N), PRIVATE :: sn_icb ! information about the calving file to be read 40 41 PUBLIC icb_init ! routine called in nemogcm.F90 module 42 PUBLIC icb_gen ! routine called in icbclv.F90 module 43 PRIVATE icb_nam 44 PRIVATE icb_alloc 45 35 PUBLIC icb_init ! routine called in nemogcm.F90 module 36 PUBLIC icb_gen ! routine called in icbclv.F90 module 37 38 CHARACTER(len=100) :: cn_dir = './' ! Root directory for location of icb files 39 TYPE(FLD_N) :: sn_icb ! information about the calving file to be read 40 41 !!---------------------------------------------------------------------- 42 !! NEMO/OPA 3.3 , NEMO Consortium (2011) 43 !! $Id:$ 44 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 45 !!---------------------------------------------------------------------- 46 46 CONTAINS 47 47 … … 52 52 ! 53 53 INTEGER :: ill 54 !55 ! open ascii output file or files for iceberg status information56 ! note that we choose to do this on all processors since we cannot predict where57 ! icebergs will be ahead of time58 !59 CALL ctl_opn( numicb, 'icebergs.stat', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', &60 -1, numout, lwp, narea )61 54 62 55 icb_alloc = 0 … … 129 122 !! ** Method : - blah blah 130 123 !!---------------------------------------------------------------------- 131 REAL(wp) , INTENT(in) :: pdt! iceberg time-step (rdt*nn_fsbc)132 INTEGER , INTENT(in) :: kt! time step number133 ! 134 INTEGER :: ji, jj, jn ! loop counters135 INTEGER :: i1, i2, i3 ! dummyintegers136 INTEGER :: ii, inum, ivar137 INTEGER :: istat1, istat2, istat3138 CHARACTER(len=300) :: cl_sdist139 !!---------------------------------------------------------------------- 140 124 REAL(wp), INTENT(in) :: pdt ! iceberg time-step (rdt*nn_fsbc) 125 INTEGER , INTENT(in) :: kt ! time step number 126 ! 127 INTEGER :: ji, jj, jn ! dummy loop indices 128 INTEGER :: i1, i2, i3 ! local integers 129 INTEGER :: ii, inum, ivar ! - - 130 INTEGER :: istat1, istat2, istat3 ! - - 131 CHARACTER(len=300) :: cl_sdist ! local character 132 !!---------------------------------------------------------------------- 133 ! 141 134 CALL icb_nam ! Read and print namelist parameters 142 IF( .NOT. ln_icebergs ) RETURN 135 ! 136 IF( .NOT. ln_icebergs ) RETURN 143 137 144 138 ! ! allocate gridded fields 145 139 IF( icb_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'icb_alloc : unable to allocate arrays' ) 146 140 141 ! ! open ascii output file or files for iceberg status information 142 ! ! note that we choose to do this on all processors since we cannot 143 ! ! predict where icebergs will be ahead of time 144 CALL ctl_opn( numicb, 'icebergs.stat', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea ) 145 147 146 ! set parameters (mostly from namelist) 148 147 ! 149 berg_dt 148 berg_dt = pdt 150 149 first_width (:) = SQRT( rn_initial_mass(:) / ( rn_LoW_ratio * rn_rho_bergs * rn_initial_thickness(:) ) ) 151 150 first_length(:) = rn_LoW_ratio * first_width(:) … … 158 157 berg_grid%stored_ice (:,:,:) = 0._wp 159 158 berg_grid%tmp (:,:) = 0._wp 160 src_calving (:,:) = 0._wp 161 src_calving_hflx (:,:) = 0._wp 162 163 ! domain for icebergs 164 IF( lk_mpp .AND. jpni == 1 ) & 165 CALL ctl_stop('icbinit: having ONE processor in x currently does not work') 166 159 src_calving (:,:) = 0._wp 160 src_calving_hflx (:,:) = 0._wp 161 162 ! ! domain for icebergs 163 IF( lk_mpp .AND. jpni == 1 ) CALL ctl_stop( 'icbinit: having ONE processor in x currently does not work' ) 164 ! NB: the issue here is simply that cyclic east-west boundary condition have not been coded in mpp case 167 165 ! for the north fold we work out which points communicate by asking 168 166 ! lbc_lnk to pass processor number (valid even in single processor case) … … 171 169 ! pack i and j together using a scaling of a power of 10 172 170 nicbpack = 10000 173 IF( jpiglo >= nicbpack ) CALL ctl_stop('icbini: processor index packing failure')171 IF( jpiglo >= nicbpack ) CALL ctl_stop( 'icbini: processor index packing failure' ) 174 172 nicbfldproc(:) = -1 175 173 … … 177 175 DO ji = 1, jpi 178 176 src_calving_hflx(ji,jj) = narea 179 src_calving(ji,jj) = nicbpack*(njmpp+jj-1) + nimpp+ji-1 180 ENDDO 181 ENDDO 177 src_calving (ji,jj) = nicbpack * (njmpp+jj-1) + nimpp+ji-1 !!gm use here mjg(jj) and mig(ji)... 178 !!gm src_calving (ji,jj) = nicbpack * mjg(jj) + mig(ji) 179 END DO 180 END DO 182 181 CALL lbc_lnk( src_calving_hflx, 'T', 1._wp ) 183 182 CALL lbc_lnk( src_calving , 'T', 1._wp ) … … 189 188 nicbdi = -1 190 189 nicbei = -1 191 DO ji = 1, jpi190 DO ji = 1, jpi 192 191 i3 = INT( src_calving(ji,jj) ) 193 192 i2 = INT( i3/nicbpack ) … … 195 194 i3 = INT( src_calving_hflx(ji,jj) ) 196 195 IF( i1 == nimpp+ji-1 .AND. i3 == narea ) THEN 197 IF( nicbdi < 0 ) THEN 198 nicbdi = ji 199 ELSE 200 nicbei = ji 196 IF( nicbdi < 0 ) THEN ; nicbdi = ji 197 ELSE ; nicbei = ji 201 198 ENDIF 202 199 ENDIF 203 200 END DO 204 201 ! 205 202 ! repeat for j direction 206 203 ji = jpi/2 207 204 nicbdj = -1 208 205 nicbej = -1 209 DO jj = 1, jpj206 DO jj = 1, jpj 210 207 i3 = INT( src_calving(ji,jj) ) 211 208 i2 = INT( i3/nicbpack ) … … 213 210 i3 = INT( src_calving_hflx(ji,jj) ) 214 211 IF( i2 == njmpp+jj-1 .AND. i3 == narea ) THEN 215 IF( nicbdj < 0 ) THEN 216 nicbdj = jj 217 ELSE 218 nicbej = jj 212 IF( nicbdj < 0 ) THEN ; nicbdj = jj 213 ELSE ; nicbej = jj 219 214 ENDIF 220 215 ENDIF 221 216 END DO 222 217 ! 223 218 ! special for east-west boundary exchange we save the destination index 224 219 i1 = MAX( nicbdi-1, 1) … … 244 239 IF( nicbfldproc(jn) == -1 ) THEN 245 240 nicbfldproc(jn) = ii 246 EXIT 241 EXIT !!gm EXIT should be avoided: use DO WHILE expression instead 247 242 ENDIF 248 243 IF( nicbfldproc(jn) == ii ) EXIT 249 END DO250 END DO244 END DO 245 END DO 251 246 ENDIF 252 247 ! … … 281 276 ENDIF 282 277 283 src_calving (:,:)= 0._wp278 src_calving (:,:) = 0._wp 284 279 src_calving_hflx(:,:) = 0._wp 285 280 … … 317 312 318 313 IF( .NOT.ln_rstart ) THEN 319 IF( nn_test_icebergs > 0 ) CALL icb_gen()314 IF( nn_test_icebergs > 0 ) CALL icb_gen() 320 315 ELSE 321 316 IF( nn_test_icebergs > 0 ) THEN 322 317 CALL icb_gen() 323 318 ELSE 324 !325 319 CALL icebergs_read_restart() 326 320 l_restarted_bergs = .TRUE. … … 337 331 338 332 SUBROUTINE icb_gen() 339 340 ! Local variables 333 !!---------------------------------------------------------------------- 334 !! *** ROUTINE icb_gen *** 335 !! 336 !! ** Purpose : iceberg generation 337 !! 338 !! ** Method : - blah blah 339 !!---------------------------------------------------------------------- 341 340 INTEGER :: ji, jj, ibergs 342 341 TYPE(iceberg) :: localberg ! NOT a pointer but an actual local variable … … 344 343 INTEGER :: iyr, imon, iday, ihr, imin, isec 345 344 INTEGER :: iberg 345 !!---------------------------------------------------------------------- 346 346 347 347 ! For convenience … … 350 350 ! call get_date(Time, iyr, imon, iday, ihr, imin, isec) 351 351 ! Convert nemo time variables from dom_oce into local versions 352 iyr = nyear352 iyr = nyear 353 353 imon = nmonth 354 354 iday = nday … … 361 361 ! use tmask here because tmask_i has been doctored on one side of the north fold line 362 362 363 DO jj = nicbdj, nicbej364 DO ji = nicbdi, nicbei365 IF (tmask(ji,jj,1) .GT. 0._wp .AND.&366 gphit(ji,jj) .GT. rn_test_box(3) .AND. gphit(ji,jj) .LT. rn_test_box(4) .AND.&367 glamt(ji,jj) .GT. rn_test_box(1) .AND. glamt(ji,jj) .LT. rn_test_box(2)) THEN363 DO jj = nicbdj, nicbej 364 DO ji = nicbdi, nicbei 365 IF( tmask(ji,jj,1) > 0._wp .AND. & 366 rn_test_box(1) < glamt(ji,jj) .AND. glamt(ji,jj) < rn_test_box(2) .AND. & 367 rn_test_box(3) < gphit(ji,jj) .AND. gphit(ji,jj) < rn_test_box(4) ) THEN 368 368 localberg%mass_scaling = rn_mass_scaling(iberg) 369 369 localpt%xi = REAL( nimpp+ji-1, wp ) … … 371 371 localpt%lon = bilin(glamt, localpt%xi, localpt%yj, 'T', 0, 0 ) 372 372 localpt%lat = bilin(gphit, localpt%xi, localpt%yj, 'T', 0, 0 ) 373 localpt%mass = rn_initial_mass(iberg)373 localpt%mass = rn_initial_mass (iberg) 374 374 localpt%thickness = rn_initial_thickness(iberg) 375 localpt%width = first_width(iberg)375 localpt%width = first_width (iberg) 376 376 localpt%length = first_length(iberg) 377 377 localpt%year = iyr … … 385 385 call add_new_berg_to_list(localberg, localpt) 386 386 ENDIF 387 END DO388 END DO389 387 END DO 388 END DO 389 ! 390 390 ibergs = count_bergs() 391 391 IF( lk_mpp ) CALL mpp_sum(ibergs) 392 392 WRITE(numicb,'(a,i6,a)') 'diamonds, icb_gen: ',ibergs,' were generated' 393 393 ! 394 394 END SUBROUTINE icb_gen 395 395 … … 398 398 !! *** ROUTINE icb_nam *** 399 399 !! 400 !! ** Purpose : read domaine namelistsand print the variables.400 !! ** Purpose : read iceberg namelist and print the variables. 401 401 !! 402 402 !! ** input : - namberg namelist 403 403 !!---------------------------------------------------------------------- 404 405 NAMELIST/namberg/ ln_icebergs, ln_bergdia, nn_sample_rate, rn_initial_mass, &406 & rn_distribution, rn_mass_scaling, rn_initial_thickness, nn_verbose_write, &407 & rn_rho_bergs, rn_LoW_ratio, nn_verbose_level, ln_operator_splitting,&408 & rn_ bits_erosion_fraction, rn_sicn_shift, ln_passive_mode,&409 & ln_time_average_weight, nn_test_icebergs, rn_test_box, rn_speed_limit,&410 & cn_dir, sn_icb411 INTEGER :: jn ! loop counter412 REAL(wp) :: zfact404 INTEGER :: jn ! dummy loop indices 405 REAL(wp) :: zfact ! local scalar 406 ! 407 NAMELIST/namberg/ ln_icebergs , ln_bergdia , nn_sample_rate , rn_initial_mass , & 408 & rn_distribution, rn_mass_scaling, rn_initial_thickness, nn_verbose_write , & 409 & rn_rho_bergs , rn_LoW_ratio , nn_verbose_level , ln_operator_splitting, & 410 & rn_bits_erosion_fraction , rn_sicn_shift , ln_passive_mode , & 411 & ln_time_average_weight , nn_test_icebergs , rn_test_box , & 412 & rn_speed_limit , cn_dir, sn_icb 413 413 !!---------------------------------------------------------------------- 414 414 … … 420 420 REWIND( numnam ) ! Namelist namrun : iceberg parameters 421 421 READ ( numnam, namberg ) 422 423 IF( .NOT. ln_icebergs .AND. lwp ) THEN 424 WRITE(numout,*) 425 WRITE(numout,*) 'icb_nam : iceberg initialization through namberg namelist read' 426 WRITE(numout,*) '~~~~~~~~~~~~~ ' 427 WRITE(numout,*) 'NO icebergs used' 422 423 !!gm bug here in mpp : the return instruction will be done only on the 1st processor (lwp=True) 424 ! IF( .NOT. ln_icebergs .AND. lwp ) THEN 425 ! WRITE(numout,*) 426 ! WRITE(numout,*) 'icb_nam : iceberg initialization through namberg namelist read' 427 ! WRITE(numout,*) '~~~~~~~~ ' 428 ! WRITE(numout,*) 'NO icebergs used' 429 ! RETURN 430 ! ENDIF 431 !!gm solution: 432 433 IF( .NOT. ln_icebergs ) THEN ! no icebergs 434 IF(lwp) THEN 435 WRITE(numout,*) 436 WRITE(numout,*) 'icb_nam : ln_icebergs = F , NO icebergs used' 437 WRITE(numout,*) '~~~~~~~~ ' 438 ENDIF 428 439 RETURN 429 440 ENDIF 441 !!gm end 430 442 431 443 IF( nn_test_icebergs > nclasses ) THEN 432 IF( lwp ) WRITE(numout,*) 'Resetting nn_test_icebergs to ',nclasses444 IF(lwp) WRITE(numout,*) 'Resetting nn_test_icebergs to ', nclasses 433 445 nn_test_icebergs = nclasses 434 446 ENDIF … … 437 449 IF( zfact < 1._wp ) THEN 438 450 IF( zfact <= 0._wp ) THEN 439 CALL ctl_stop( 'icb_init: sum of berg distribution equal to zero')451 CALL ctl_stop( 'icb_init: sum of berg distribution equal to zero' ) 440 452 ELSE 441 453 rn_distribution(:) = rn_distribution(:) / zfact 442 CALL ctl_warn( 'icb_init: sum of berg input distribution not equal to one and so RESCALED')454 CALL ctl_warn( 'icb_init: sum of berg input distribution not equal to one and so RESCALED' ) 443 455 ENDIF 444 456 ENDIF … … 447 459 WRITE(numout,*) 448 460 WRITE(numout,*) 'icb_nam : iceberg initialization through namberg namelist read' 449 WRITE(numout,*) '~~~~~~~~ ~~~~~'461 WRITE(numout,*) '~~~~~~~~ ' 450 462 WRITE(numout,*) ' Calculate budgets ln_bergdia = ', ln_bergdia 451 463 WRITE(numout,*) ' Period between sampling of position for trajectory storage nn_sample_rate = ', nn_sample_rate … … 455 467 ENDDO 456 468 WRITE(numout,*) ' Fraction of calving to apply to this class (non-dim) rn_distribution =' 457 DO jn =1,nclasses469 DO jn = 1, nclasses 458 470 WRITE(numout,'(a,f10.2)') ' ',rn_distribution(jn) 459 END DO471 END DO 460 472 WRITE(numout,*) ' Ratio between effective and real iceberg mass (non-dim) rn_mass_scaling = ' 461 DO jn =1,nclasses473 DO jn = 1, nclasses 462 474 WRITE(numout,'(a,f10.2)') ' ',rn_mass_scaling(jn) 463 END DO464 WRITE(numout,*) ' Total thickness of newly calved bergs (m) rn_initial_thickness = '465 DO jn =1,nclasses475 END DO 476 WRITE(numout,*) ' Total thickness of newly calved bergs (m) rn_initial_thickness = ' 477 DO jn = 1, nclasses 466 478 WRITE(numout,'(a,f10.2)') ' ',rn_initial_thickness(jn) 467 END DO479 END DO 468 480 WRITE(numout,*) ' Timesteps between verbose messages nn_verbose_write = ', nn_verbose_write 469 481 470 WRITE(numout,*) ' Density of icebergs rn_rho_bergs 471 WRITE(numout,*) ' Initial ratio L/W for newly calved icebergs rn_LoW_ratio 472 WRITE(numout,*) ' Turn on more verbose output level= ', nn_verbose_level482 WRITE(numout,*) ' Density of icebergs rn_rho_bergs = ', rn_rho_bergs 483 WRITE(numout,*) ' Initial ratio L/W for newly calved icebergs rn_LoW_ratio = ', rn_LoW_ratio 484 WRITE(numout,*) ' Turn on more verbose output level = ', nn_verbose_level 473 485 WRITE(numout,*) ' Use first order operator splitting for thermodynamics ', & 474 486 & 'use_operator_splitting = ', ln_operator_splitting … … 477 489 478 490 WRITE(numout,*) ' Shift of sea-ice concentration in erosion flux modulation ', & 479 & '(0<sicn_shift<1) sicn_shift= ', rn_sicn_shift491 & '(0<sicn_shift<1) rn_sicn_shift = ', rn_sicn_shift 480 492 WRITE(numout,*) ' Do not add freshwater flux from icebergs to ocean ', & 481 493 & ' passive_mode = ', ln_passive_mode 482 WRITE(numout,*) ' Time average the weight on the ocean time_average_weight = ', ln_time_average_weight494 WRITE(numout,*) ' Time average the weight on the ocean time_average_weight = ', ln_time_average_weight 483 495 WRITE(numout,*) ' Create icebergs in absence of a restart file nn_test_icebergs = ', nn_test_icebergs 484 WRITE(numout,*) ' in lon/lat box 485 WRITE(numout,*) ' CFL speed limit for a berg speed_limit = ', rn_speed_limit496 WRITE(numout,*) ' in lon/lat box = ', rn_test_box 497 WRITE(numout,*) ' CFL speed limit for a berg speed_limit = ', rn_speed_limit 486 498 WRITE(numout,*) ' Writing Iceberg status information to icebergs.stat file ' 487 499 ENDIF -
branches/2012/dev_r3337_NOCS10_ICB/NEMOGCM/NEMO/OPA_SRC/ICB/icblbc.F90
r3360 r3370 75 75 #endif 76 76 77 !!------------------------------------------------------------------------- 77 !!---------------------------------------------------------------------- 78 !! NEMO/OPA 3.3 , NEMO Consortium (2011) 79 !! $Id:$ 80 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 81 !!---------------------------------------------------------------------- 78 82 CONTAINS 79 83 … … 85 89 !! including north-fold 86 90 !!---------------------------------------------------------------------- 87 ! Local variables88 TYPE( iceberg), POINTER :: this89 TYPE(point) , POINTER :: pt90 INTEGER :: iine91 TYPE(iceberg), POINTER :: this 92 TYPE(point) , POINTER :: pt 93 INTEGER :: iine 94 !!---------------------------------------------------------------------- 91 95 92 96 !! periodic east/west boundaries … … 99 103 pt => this%current_point 100 104 iine = INT( pt%xi + 0.5 ) 101 IF( iine .GT.nimpp+nicbei-1 ) THEN105 IF( iine > nimpp+nicbei-1 ) THEN 102 106 pt%xi = ricb_right + MOD(pt%xi, 1._wp ) - 1._wp 103 ELSE IF( iine .LT.nimpp+nicbdi-1 ) THEN107 ELSE IF( iine < nimpp+nicbdi-1 ) THEN 104 108 pt%xi = ricb_left + MOD(pt%xi, 1._wp ) 105 109 ENDIF 106 110 this => this%next 107 END DO108 111 END DO 112 ! 109 113 ENDIF 110 114 111 115 !! north/south boundaries 112 116 !! ====================== 113 114 !! south symmetric 115 116 IF( nperio == 2 ) THEN 117 CALL ctl_stop(' south symmetric condition not implemented for icebergs') 118 ENDIF 119 120 !! north fold 121 122 IF( nperio == 3 .OR. nperio == 4 .OR. nperio == 5 .OR. nperio == 6 ) CALL lbc_nfld_bergs() 123 117 ! south symmetric 118 IF( nperio == 2 ) CALL ctl_stop(' south symmetric condition not implemented for icebergs') 119 ! north fold 120 IF( nperio == 3 .OR. nperio == 4 .OR. nperio == 5 .OR. nperio == 6 ) CALL lbc_nfld_bergs() 121 ! 124 122 END SUBROUTINE lbc_send_bergs 125 123 126 !!-------------------------------------------------------------------------127 124 128 125 SUBROUTINE lbc_nfld_bergs() … … 132 129 !! ** Purpose : single processor north fold exchange 133 130 !!---------------------------------------------------------------------- 134 ! 135 ! Local variables 136 TYPE(iceberg), POINTER :: this 137 TYPE(point) , POINTER :: pt 138 INTEGER :: iine, ijne, ipts 139 INTEGER :: iiglo, ijglo 140 131 TYPE(iceberg), POINTER :: this 132 TYPE(point) , POINTER :: pt 133 INTEGER :: iine, ijne, ipts 134 INTEGER :: iiglo, ijglo 135 !!---------------------------------------------------------------------- 136 ! 141 137 this => first_berg 142 138 DO WHILE( ASSOCIATED(this) ) … … 158 154 ENDIF 159 155 this => this%next 160 END DO161 156 END DO 157 ! 162 158 END SUBROUTINE lbc_nfld_bergs 163 159 164 !!-------------------------------------------------------------------------165 166 160 #if defined key_mpp_mpi 161 !!---------------------------------------------------------------------- 162 !! 'key_mpp_mpi' MPI massively parallel processing library 163 !!---------------------------------------------------------------------- 167 164 168 165 SUBROUTINE mpp_send_bergs() … … 177 174 !! then buffer is sent if necessary 178 175 !!---------------------------------------------------------------------- 179 180 ! Local variables 181 TYPE(iceberg) , POINTER :: tmpberg, this 182 TYPE(point) , POINTER :: pt 183 INTEGER :: ibergs_to_send_e, ibergs_to_send_w 184 INTEGER :: ibergs_to_send_n, ibergs_to_send_s 185 INTEGER :: ibergs_rcvd_from_e, ibergs_rcvd_from_w 186 INTEGER :: ibergs_rcvd_from_n, ibergs_rcvd_from_s 187 INTEGER :: i, ibergs_start, ibergs_end 188 INTEGER :: iine, ijne 189 INTEGER :: ipe_N, ipe_S, ipe_W, ipe_E 190 REAL(wp), DIMENSION(2) :: zewbergs, zwebergs, znsbergs, zsnbergs 191 INTEGER :: iml_req1, iml_req2, iml_req3, iml_req4 192 INTEGER :: iml_req5, iml_req6, iml_req7, iml_req8, iml_err 193 INTEGER, DIMENSION(MPI_STATUS_SIZE) :: iml_stat 176 TYPE(iceberg) , POINTER :: tmpberg, this 177 TYPE(point) , POINTER :: pt 178 INTEGER :: ibergs_to_send_e, ibergs_to_send_w 179 INTEGER :: ibergs_to_send_n, ibergs_to_send_s 180 INTEGER :: ibergs_rcvd_from_e, ibergs_rcvd_from_w 181 INTEGER :: ibergs_rcvd_from_n, ibergs_rcvd_from_s 182 INTEGER :: i, ibergs_start, ibergs_end 183 INTEGER :: iine, ijne 184 INTEGER :: ipe_N, ipe_S, ipe_W, ipe_E 185 REAL(wp), DIMENSION(2) :: zewbergs, zwebergs, znsbergs, zsnbergs 186 INTEGER :: iml_req1, iml_req2, iml_req3, iml_req4 187 INTEGER :: iml_req5, iml_req6, iml_req7, iml_req8, iml_err 188 INTEGER, DIMENSION(MPI_STATUS_SIZE) :: iml_stat 194 189 195 190 ! set up indices of neighbouring processors … … 237 232 ENDIF 238 233 239 ibergs_to_send_e = 0240 ibergs_to_send_w = 0241 ibergs_to_send_n = 0242 ibergs_to_send_s = 0234 ibergs_to_send_e = 0 235 ibergs_to_send_w = 0 236 ibergs_to_send_n = 0 237 ibergs_to_send_s = 0 243 238 ibergs_rcvd_from_e = 0 244 239 ibergs_rcvd_from_w = 0 … … 246 241 ibergs_rcvd_from_s = 0 247 242 248 ! Find number of bergs that headed east/west 249 IF( ASSOCIATED(first_berg) ) THEN 243 IF( ASSOCIATED(first_berg) ) THEN ! Find number of bergs that headed east/west 250 244 this => first_berg 251 245 DO WHILE (ASSOCIATED(this)) 252 246 pt => this%current_point 253 247 iine = INT( pt%xi + 0.5 ) 254 IF( ipe_E >= 0 .AND. iine .GT.nimpp+nicbei-1 ) THEN248 IF( ipe_E >= 0 .AND. iine > nimpp+nicbei-1 ) THEN 255 249 tmpberg => this 256 250 this => this%next … … 265 259 CALL pack_berg_into_buffer( tmpberg, obuffer_e, ibergs_to_send_e) 266 260 CALL delete_iceberg_from_list(first_berg, tmpberg) 267 ELSE IF( ipe_W >= 0 .AND. iine .LT.nimpp+nicbdi-1 ) THEN261 ELSE IF( ipe_W >= 0 .AND. iine < nimpp+nicbdi-1 ) THEN 268 262 tmpberg => this 269 263 this => this%next … … 281 275 this => this%next 282 276 ENDIF 283 ENDDO 284 ENDIF 285 if( nn_verbose_level >= 3) then 286 write(numicb,*) 'bergstep ',nktberg,' send ew: ', ibergs_to_send_e, ibergs_to_send_w 287 call flush(numicb) 288 endif 289 290 ! send bergs east and receive bergs from west (ie ones that were sent east) 291 ! and vice versa 277 END DO 278 ENDIF 279 IF( nn_verbose_level >= 3) THEN 280 WRITE(numicb,*) 'bergstep ',nktberg,' send ew: ', ibergs_to_send_e, ibergs_to_send_w 281 CALL flush(numicb) 282 ENDIF 283 284 ! send bergs east and receive bergs from west (ie ones that were sent east) and vice versa 292 285 293 286 ! pattern here is copied from lib_mpp code … … 318 311 ibergs_rcvd_from_w = INT( zwebergs(2) ) 319 312 END SELECT 320 if( nn_verbose_level >= 3) then321 write(numicb,*) 'bergstep ',nktberg,' recv ew: ', ibergs_rcvd_from_w, ibergs_rcvd_from_e322 callflush(numicb)323 endif313 IF( nn_verbose_level >= 3) THEN 314 WRITE(numicb,*) 'bergstep ',nktberg,' recv ew: ', ibergs_rcvd_from_w, ibergs_rcvd_from_e 315 CALL flush(numicb) 316 ENDIF 324 317 325 318 SELECT CASE ( nbondi ) … … 357 350 ENDIF 358 351 CALL unpack_berg_from_buffer(first_berg, ibuffer_e, i) 359 END DO352 END DO 360 353 DO i = 1, ibergs_rcvd_from_w 361 354 IF( nn_verbose_level >= 4 ) THEN … … 378 371 ENDIF 379 372 CALL unpack_berg_from_buffer(first_berg, ibuffer_w, i) 380 END DO373 END DO 381 374 END SELECT 382 375 … … 414 407 this => this%next 415 408 ENDIF 416 END DO409 END DO 417 410 ENDIF 418 411 if( nn_verbose_level >= 3) then … … 468 461 ENDIF 469 462 CALL unpack_berg_from_buffer(first_berg, ibuffer_n, i) 470 END DO463 END DO 471 464 CASE( 0 ) 472 465 IF( ibergs_to_send_s > 0 ) CALL mppsend( 17, obuffer_s%data, ibergs_to_send_s*jp_buffer_width, ipe_S, iml_req2 ) … … 488 481 ENDIF 489 482 CALL unpack_berg_from_buffer(first_berg, ibuffer_n, i) 490 END DO483 END DO 491 484 DO i = 1, ibergs_rcvd_from_s 492 485 IF( nn_verbose_level >= 4 ) THEN … … 509 502 ENDIF 510 503 CALL unpack_berg_from_buffer(first_berg, ibuffer_s, i) 511 END DO504 END DO 512 505 END SELECT 513 506 … … 578 571 ENDIF ! root_pe 579 572 ENDIF ! debug 580 573 ! 581 574 CALL mppsync() 582 575 ! 583 576 END SUBROUTINE mpp_send_bergs 584 577 585 !!-------------------------------------------------------------------------586 578 587 579 SUBROUTINE mpp_nfld_bergs() 588 589 ! Local variables 580 !!---------------------------------------------------------------------- 581 !! *** SUBROUTINE mpp_nfld_bergs *** 582 !! 583 !! ** Purpose : north fold treatment in multi processor exchange 584 !! 585 !! ** Method : 586 !!---------------------------------------------------------------------- 590 587 TYPE(iceberg) , POINTER :: tmpberg, this 591 588 TYPE(point) , POINTER :: pt … … 649 646 ENDIF 650 647 this => this%next 651 END DO648 END DO 652 649 ENDIF 653 650 if( nn_verbose_level >= 3) then … … 681 678 ENDIF 682 679 CALL unpack_berg_from_buffer(first_berg, ibuffer_f, jk ) 683 END DO684 685 END DO686 680 END DO 681 ! 682 END DO 683 ! 687 684 END SUBROUTINE mpp_nfld_bergs 688 685 689 !!-------------------------------------------------------------------------690 686 691 687 SUBROUTINE dealloc_buffers() 692 688 !!---------------------------------------------------------------------- 689 !!---------------------------------------------------------------------- 693 690 CALL dealloc_buffer( obuffer_n ) 694 691 CALL dealloc_buffer( obuffer_s ) … … 702 699 END SUBROUTINE dealloc_buffers 703 700 704 !!------------------------------------------------------------------------- 705 706 SUBROUTINE dealloc_buffer(pbuff)707 ! Arguments701 702 SUBROUTINE dealloc_buffer( pbuff ) 703 !!---------------------------------------------------------------------- 704 !!---------------------------------------------------------------------- 708 705 TYPE(buffer), POINTER :: pbuff 709 706 !!---------------------------------------------------------------------- 710 707 IF( ASSOCIATED(pbuff) ) THEN 711 708 IF( ASSOCIATED(pbuff%data)) DEALLOCATE(pbuff%data) … … 714 711 END SUBROUTINE dealloc_buffer 715 712 716 !!------------------------------------------------------------------------- 717 718 SUBROUTINE pack_berg_into_buffer(berg, pbuff, kb) 719 ! Arguments 720 TYPE(iceberg), POINTER :: berg 721 TYPE(buffer) , POINTER :: pbuff 722 INTEGER , INTENT(in) :: kb 723 ! Local 724 INTEGER :: k 725 713 714 SUBROUTINE pack_berg_into_buffer( berg, pbuff, kb ) 715 !!---------------------------------------------------------------------- 716 !!---------------------------------------------------------------------- 717 TYPE(iceberg), POINTER :: berg 718 TYPE(buffer) , POINTER :: pbuff 719 INTEGER , INTENT(in) :: kb 720 ! 721 INTEGER :: k ! local integer 722 !!---------------------------------------------------------------------- 723 ! 726 724 IF( .NOT. ASSOCIATED(pbuff) ) CALL increase_buffer( pbuff, jp_delta_buf ) 727 725 IF( kb .GT. pbuff%size ) CALL increase_buffer( pbuff, jp_delta_buf ) … … 748 746 pbuff%data(15+k,kb) = REAL( berg%number(k), wp ) 749 747 END DO 750 748 ! 751 749 END SUBROUTINE pack_berg_into_buffer 752 750 753 !!-------------------------------------------------------------------------754 751 755 752 SUBROUTINE unpack_berg_from_buffer(first, pbuff, kb) 756 ! Arguments 753 !!---------------------------------------------------------------------- 754 !!---------------------------------------------------------------------- 757 755 TYPE(iceberg), POINTER :: first 758 756 TYPE(buffer) , POINTER :: pbuff 759 757 INTEGER , INTENT(in) :: kb 760 ! Local variables758 ! 761 759 TYPE(iceberg) :: currentberg 762 760 TYPE(point) :: pt 763 761 INTEGER :: ik 764 762 !!---------------------------------------------------------------------- 763 ! 765 764 pt%lon = pbuff%data( 1,kb) 766 765 pt%lat = pbuff%data( 2,kb) … … 779 778 780 779 currentberg%mass_scaling = pbuff%data(15,kb) 781 DO ik =1,nkounts780 DO ik = 1, nkounts 782 781 currentberg%number(ik) = INT( pbuff%data(15+ik,kb) ) 783 782 END DO 784 783 ! 785 784 CALL add_new_berg_to_list(currentberg, pt ) 786 785 ! 787 786 END SUBROUTINE unpack_berg_from_buffer 788 787 789 !!-------------------------------------------------------------------------790 788 791 789 SUBROUTINE increase_buffer(old,kdelta) 792 ! Arguments 793 TYPE(buffer), POINTER :: old 794 INTEGER , INTENT(in) :: kdelta 795 ! Local variables 796 TYPE(buffer), POINTER :: new 797 INTEGER :: inew_size 798 799 IF( .NOT. ASSOCIATED(old) ) THEN 800 inew_size = kdelta 801 ELSE 802 inew_size = old%size + kdelta 790 !!---------------------------------------------------------------------- 791 TYPE(buffer), POINTER :: old 792 INTEGER , INTENT(in) :: kdelta 793 ! 794 TYPE(buffer), POINTER :: new 795 INTEGER :: inew_size 796 !!---------------------------------------------------------------------- 797 ! 798 IF( .NOT. ASSOCIATED(old) ) THEN ; inew_size = kdelta 799 ELSE ; inew_size = old%size + kdelta 803 800 ENDIF 804 801 ALLOCATE( new ) … … 811 808 ENDIF 812 809 old => new 813 810 ! 814 811 END SUBROUTINE increase_buffer 815 812 816 !!-------------------------------------------------------------------------817 813 818 814 SUBROUTINE increase_ibuffer(old,kdelta) 819 ! Arguments 815 !!---------------------------------------------------------------------- 816 !!---------------------------------------------------------------------- 820 817 TYPE(buffer), POINTER :: old 821 818 INTEGER , INTENT(in) :: kdelta 822 ! Local variables819 ! 823 820 TYPE(buffer), POINTER :: new 824 821 INTEGER :: inew_size, iold_size 822 !!---------------------------------------------------------------------- 825 823 826 824 IF( .NOT. ASSOCIATED(old) ) THEN … … 848 846 !WRITE( numicb,*) 'increase_ibuffer',narea,' increased to',inew_size 849 847 ENDIF 850 848 ! 851 849 END SUBROUTINE increase_ibuffer 852 850 853 !!-------------------------------------------------------------------------854 855 851 #else 856 852 !!---------------------------------------------------------------------- 853 !! Default case: Dummy module share memory computing 854 !!---------------------------------------------------------------------- 857 855 SUBROUTINE mpp_send_bergs() 858 856 WRITE(numout,*) 'mpp_send_bergs: You should not have seen this message!!' … … 865 863 #endif 866 864 865 !!====================================================================== 867 866 END MODULE icblbc -
branches/2012/dev_r3337_NOCS10_ICB/NEMOGCM/NEMO/OPA_SRC/ICB/icbrst.F90
r3361 r3370 14 14 !!---------------------------------------------------------------------- 15 15 !!---------------------------------------------------------------------- 16 !! icebergs_read_restart : initialise 17 !! icebergs_write_restart : generate test icebergs 16 !! icebergs_read_restart : initialise !!gm suggested name : icebergs_rst_read or better icb_rst_read 17 !! icebergs_write_restart : generate test icebergs !!gm icebergs_rst_write or better icb_rst_write 18 18 !!---------------------------------------------------------------------- 19 19 USE par_oce ! NEMO parameters … … 28 28 PRIVATE 29 29 30 REAL(wp), DIMENSION(:,:,:), PRIVATE, POINTER :: griddata => NULL() ! need 2d array to read in with 31 INTEGER, DIMENSION(3), PRIVATE :: nstrt3, nlngth3 32 INTEGER, PRIVATE :: nlonid, nlatid, nxid, nyid, nuvelid, nvvelid 33 INTEGER, PRIVATE :: nmassid, nthicknessid, nwidthid, nlengthid 34 INTEGER, PRIVATE :: nyearid, ndayid 35 INTEGER, PRIVATE :: nscaling_id, nmass_of_bits_id, nheat_density_id, numberid 36 INTEGER, PRIVATE :: nsiceid, nsheatid, ncalvid, ncalvhid, nkountid 37 INTEGER, PRIVATE :: nret, ncid, nc_dim 38 39 PUBLIC icebergs_read_restart ! routine called in icbini.F90 module 40 PUBLIC icebergs_write_restart ! routine called in icbrun.F90 module 41 30 PUBLIC icebergs_read_restart ! routine called in icbini.F90 module 31 PUBLIC icebergs_write_restart ! routine called in icbrun.F90 module 32 33 INTEGER :: nlonid, nlatid, nxid, nyid, nuvelid, nvvelid 34 INTEGER :: nmassid, nthicknessid, nwidthid, nlengthid 35 INTEGER :: nyearid, ndayid 36 INTEGER :: nscaling_id, nmass_of_bits_id, nheat_density_id, numberid 37 INTEGER :: nsiceid, nsheatid, ncalvid, ncalvhid, nkountid 38 INTEGER :: nret, ncid, nc_dim 39 40 INTEGER, DIMENSION(3) :: nstrt3, nlngth3 41 REAL(wp), DIMENSION(:,:,:), POINTER :: griddata => NULL() ! need 2d array to read in with 42 43 !!---------------------------------------------------------------------- 44 !! NEMO/OPA 3.3 , NEMO Consortium (2011) 45 !! $Id:$ 46 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 47 !!---------------------------------------------------------------------- 42 48 CONTAINS 43 49 44 50 SUBROUTINE icebergs_read_restart() 45 ! sga - read a restart file 46 ! for this version, we just read back in the restart for this processor 47 ! so we cannot change the processor layout currently with iceberg code 48 49 ! Local variables 50 INTEGER :: idim, ivar, iatt 51 INTEGER :: jn, iunlim_dim, ibergs_in_file 52 INTEGER :: iclass 53 INTEGER, DIMENSION(1) :: istrt, ilngth, idata 54 INTEGER, DIMENSION(2) :: istrt2, ilngth2 55 INTEGER, DIMENSION(nkounts) :: idata2 56 REAL(wp), DIMENSION(1) :: zdata ! need 1d array to read in with 51 !!---------------------------------------------------------------------- 52 !! *** SUBROUTINE icebergs_read_restart *** 53 !! 54 !! ** Purpose : read a iceberg restart file 55 !! NB: for this version, we just read back in the restart for this processor 56 !! so we cannot change the processor layout currently with iceberg code 57 !!---------------------------------------------------------------------- 58 INTEGER :: idim, ivar, iatt 59 INTEGER :: jn, iunlim_dim, ibergs_in_file 60 INTEGER :: iclass 61 INTEGER, DIMENSION(1) :: istrt, ilngth, idata 62 INTEGER, DIMENSION(2) :: istrt2, ilngth2 63 INTEGER, DIMENSION(nkounts) :: idata2 64 REAL(wp), DIMENSION(1) :: zdata ! need 1d array to read in with 57 65 ! start and count arrays 58 LOGICAL :: ll_found_restart, ll_multiPErestart=.FALSE. 59 CHARACTER(len=80) :: cl_filename 60 CHARACTER(len=NF90_MAX_NAME) :: cl_dname 61 TYPE(iceberg) :: localberg ! NOT a pointer but an actual local variable 62 TYPE(point) :: localpt ! NOT a pointer but an actual local variable 66 LOGICAL :: ll_found_restart, ll_multiPErestart=.FALSE. 67 CHARACTER(len=80) :: cl_filename 68 CHARACTER(len=NF90_MAX_NAME) :: cl_dname 69 TYPE(iceberg) :: localberg ! NOT a pointer but an actual local variable 70 TYPE(point) :: localpt ! NOT a pointer but an actual local variable 71 !!---------------------------------------------------------------------- 63 72 64 73 IF( .NOT. ASSOCIATED(griddata) ) ALLOCATE( griddata(jpi,jpj,1) ) … … 168 177 ! 169 178 CALL add_new_berg_to_list( localberg, localpt ) 170 END DO171 179 END DO 180 ! 172 181 ENDIF 173 182 … … 190 199 nlngth3(3) = 1 191 200 192 DO jn=1,iclass 193 201 DO jn = 1, iclass 194 202 nstrt3(3) = jn 195 196 nret = NF90_GET_VAR( ncid, nsiceid , griddata, nstrt3, nlngth3 ) 203 nret = NF90_GET_VAR( ncid, nsiceid , griddata, nstrt3, nlngth3 ) 197 204 berg_grid%stored_ice(:,:,jn) = griddata(:,:,1) 198 199 ENDDO 200 201 nret = NF90_GET_VAR( ncid, ncalvid , src_calving(:,:) ) 202 nret = NF90_GET_VAR( ncid, ncalvhid, src_calving_hflx(:,:) ) 205 END DO 206 207 nret = NF90_GET_VAR( ncid, ncalvid , src_calving (:,:) ) 208 nret = NF90_GET_VAR( ncid, ncalvhid, src_calving_hflx (:,:) ) 203 209 nret = NF90_GET_VAR( ncid, nsheatid, berg_grid%stored_heat(:,:) ) 204 210 nret = NF90_GET_VAR( ncid, nkountid, idata2(:) ) … … 217 223 CALL mpp_sum(jn) 218 224 ENDIF 219 IF ( lwp ) THEN 220 WRITE(numout,'(a,i5,a,i5,a)') 'icebergs, read_restart_bergs: there were',ibergs_in_file, & 221 ' bergs in the restart file and', jn,' bergs have been read' 222 ENDIF 223 224 IF ( lwp .and. nn_verbose_level >= 0) & 225 WRITE(numout,'(a)') 'icebergs, read_restart_bergs: completed' 226 225 IF(lwp) WRITE(numout,'(a,i5,a,i5,a)') 'icebergs, read_restart_bergs: there were',ibergs_in_file, & 226 & ' bergs in the restart file and', jn,' bergs have been read' 227 ! 228 IF( lwp .and. nn_verbose_level >= 0) WRITE(numout,'(a)') 'icebergs, read_restart_bergs: completed' 229 ! 227 230 END SUBROUTINE icebergs_read_restart 228 231 229 ! ##############################################################################230 232 231 233 SUBROUTINE icebergs_write_restart( kt ) 232 233 INTEGER, INTENT( in ) :: kt 234 ! Local variables 235 INTEGER :: ix_dim, iy_dim, ik_dim, in_dim 236 INTEGER :: jn 237 CHARACTER(len=80) :: cl_filename 238 TYPE(iceberg), POINTER :: this 239 TYPE(point) , POINTER :: pt 234 !!---------------------------------------------------------------------- 235 !! *** SUBROUTINE icebergs_write_restart *** 236 !! 237 !!---------------------------------------------------------------------- 238 INTEGER, INTENT( in ) :: kt 239 ! 240 INTEGER :: jn ! dummy loop index 241 INTEGER :: ix_dim, iy_dim, ik_dim, in_dim 242 CHARACTER(len=80) :: cl_filename 243 TYPE(iceberg), POINTER :: this 244 TYPE(point) , POINTER :: pt 245 !!---------------------------------------------------------------------- 240 246 241 247 IF( .NOT. ASSOCIATED(griddata) ) ALLOCATE( griddata(jpi,jpj,1) ) … … 406 412 407 413 this=>this%next 408 END DO409 414 END DO 415 ! 410 416 ENDIF ! associated(first_berg) 411 417 412 418 ! Finish up 413 419 nret = NF90_CLOSE(ncid) 414 IF (nret .ne.NF90_NOERR) CALL ctl_stop('icebergs, write_restart: nf_close failed')415 420 IF (nret /= NF90_NOERR) CALL ctl_stop('icebergs, write_restart: nf_close failed') 421 ! 416 422 END SUBROUTINE icebergs_write_restart 417 423 ! 418 424 END MODULE icbrst -
branches/2012/dev_r3337_NOCS10_ICB/NEMOGCM/NEMO/OPA_SRC/ICB/icbrun.F90
r3359 r3370 3 3 !!====================================================================== 4 4 !! *** MODULE icbrun *** 5 !! Ocean physics: initialise variables for iceberg tracking5 !! Icebergs: initialise variables for iceberg tracking 6 6 !!====================================================================== 7 7 !! History : 3.3.1 ! 2010-01 (Martin&Adcroft) Original code … … 42 42 PUBLIC icb_end ! routine called in nemogcm.F90 module 43 43 44 !!---------------------------------------------------------------------- 45 !! NEMO/OPA 3.3 , NEMO Consortium (2011) 46 !! $Id:$ 47 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 48 !!---------------------------------------------------------------------- 44 49 CONTAINS 45 50 … … 52 57 !! ** Method : - blah blah 53 58 !!---------------------------------------------------------------------- 54 INTEGER, INTENT( in ) :: kt55 ! 56 INTEGER :: iyr, imon, iday, ihr, imin, isec57 LOGICAL :: ll_sample_traj, ll_budget, ll_verbose59 INTEGER, INTENT(in) :: kt ! time step index 60 ! 61 INTEGER :: iyr, imon, iday, ihr, imin, isec ! local integers 62 LOGICAL :: ll_sample_traj, ll_budget, ll_verbose ! local logical 58 63 !!---------------------------------------------------------------------- 59 64 !! start of timestep housekeeping … … 61 66 nktberg = kt 62 67 63 ! read calving data 64 IF( nn_test_icebergs < 0 ) THEN 68 IF( nn_test_icebergs < 0 ) THEN ! read calving data 65 69 ! 66 70 CALL fld_read ( kt, 1, sf_icb ) … … 88 92 89 93 ll_verbose = .FALSE. 90 IF( nn_verbose_write .GT.0 .AND. &91 MOD(kt-1,nn_verbose_write ) == 0 ) ll_verbose = nn_verbose_level >= 094 IF( nn_verbose_write > 0 .AND. & 95 MOD(kt-1,nn_verbose_write ) == 0 ) ll_verbose = nn_verbose_level >= 0 92 96 IF( ll_verbose ) WRITE(numicb,9100) iyr, imon, iday, ihr, imin, isec, & 93 97 current_year, current_yearday … … 101 105 !! process icebergs 102 106 103 ! Accumulate ice from calving 104 CALL accumulate_calving( kt ) 105 106 ! Calve excess stored ice into icebergs 107 CALL calve_icebergs() 107 CALL accumulate_calving( kt ) ! Accumulate ice from calving 108 109 CALL calve_icebergs() ! Calve excess stored ice into icebergs 110 108 111 109 112 ! !== For each berg, evolve ==! 110 113 ! 111 IF( ASSOCIATED(first_berg) ) CALL evolve_icebergs() ! ice berg dynamics 112 113 IF( lk_mpp ) THEN 114 CALL mpp_send_bergs () ! Send bergs to other PEs 115 ELSE 116 CALL lbc_send_bergs() ! Deal with any cyclic boundaries in non-mpp case 114 IF( ASSOCIATED(first_berg) ) CALL evolve_icebergs() ! ice berg dynamics 115 116 IF( lk_mpp ) THEN ; CALL mpp_send_bergs () ! Send bergs to other PEs 117 ELSE ; CALL lbc_send_bergs() ! Deal with any cyclic boundaries in non-mpp case 117 118 ENDIF 118 119 119 IF( ASSOCIATED(first_berg) ) CALL thermodynamics ( kt ) ! Ice berg thermodynamics (melting) + rolling120 IF( ASSOCIATED(first_berg) ) CALL thermodynamics ( kt ) ! Ice berg thermodynamics (melting) + rolling 120 121 121 122 !!---------------------------------------------------------------------- … … 123 124 124 125 ll_sample_traj = .FALSE. 125 IF( nn_sample_rate .GT.0 .AND. MOD(kt-1,nn_sample_rate) == 0 ) ll_sample_traj = .TRUE.126 IF( nn_sample_rate > 0 .AND. MOD(kt-1,nn_sample_rate) == 0 ) ll_sample_traj = .TRUE. 126 127 IF( ll_sample_traj .AND. & 127 128 ASSOCIATED(first_berg) ) CALL traj_write ( kt ) ! For each berg, record trajectory … … 135 136 CALL iom_put( "berg_stored_ice" , berg_grid%stored_ice (:,:,:) ) ! 'Accumulated ice mass by class', 'kg' 136 137 137 ! write out mean budgets - not sure why this happens before they are calculated - sga !138 ! write out mean budgets 138 139 CALL icb_budget_put() 139 140 140 141 ! Dump icebergs to screen 141 if ( nn_verbose_level >= 2 ) CALL print_bergs( 'icb_stp, status', kt )142 if ( nn_verbose_level >= 2 ) CALL print_bergs( 'icb_stp, status', kt ) 142 143 143 144 ! Diagnose budgets 144 145 ll_budget = .FALSE. 145 IF( nn_verbose_write .GT.0 .AND. MOD(kt-1,nn_verbose_write) == 0 ) ll_budget = ln_bergdia146 IF( nn_verbose_write > 0 .AND. MOD(kt-1,nn_verbose_write) == 0 ) ll_budget = ln_bergdia 146 147 CALL icb_budget( ll_budget ) 147 148 148 149 IF( MOD(kt,nn_stock) == 0 ) THEN 149 150 CALL icebergs_write_restart( kt ) 150 IF( nn_sample_rate .GT. 0 )CALL traj_sync()151 IF( nn_sample_rate > 0 ) CALL traj_sync() 151 152 ENDIF 152 153 … … 154 155 END SUBROUTINE icb_stp 155 156 156 !!-------------------------------------------------------------------------157 157 158 158 SUBROUTINE icb_end( kt ) 159 160 ! Arguments 159 !!---------------------------------------------------------------------- 160 !! *** ROUTINE icb_end *** 161 !! 162 !! ** Purpose : deallocate icebergs arrays and 163 !! 164 !!---------------------------------------------------------------------- 161 165 INTEGER, INTENT( in ) :: kt 162 ! Local variables166 ! 163 167 TYPE(iceberg), POINTER :: this, next 168 !!---------------------------------------------------------------------- 164 169 165 170 ! expanded arrays for bilinear interpolation … … 195 200 CALL destroy_iceberg(this) 196 201 this=>next 197 END DO202 END DO 198 203 199 204 CALL icb_budget_end() … … 211 216 DEALLOCATE(first_length) 212 217 213 IF (lwp) WRITE(numout,'(a,i6)') 'icebergs: icb_end complete',narea218 IF(lwp) WRITE(numout,'(a,i6)') 'icebergs: icb_end complete', narea 214 219 CALL flush( numicb ) 215 220 CLOSE( numicb ) 216 221 ! 217 222 END SUBROUTINE icb_end 218 223 -
branches/2012/dev_r3337_NOCS10_ICB/NEMOGCM/NEMO/OPA_SRC/ICB/icbthm.F90
r3359 r3370 3 3 !!====================================================================== 4 4 !! *** MODULE icbthm *** 5 !! Ocean physics: thermodynamics routines for icebergs5 !! Icebergs: thermodynamics routines for icebergs 6 6 !!====================================================================== 7 7 !! History : 3.3.1 ! 2010-01 (Martin&Adcroft) Original code … … 41 41 !! ** Method : - blah blah 42 42 !!---------------------------------------------------------------------- 43 INTEGER :: kt! timestep number, just passed to print_berg44 ! 45 REAL(wp) :: zM, zT, zW, zL, zSST, zVol, zLn, zWn, zTn, znVol, zIC, zDn46 REAL(wp) :: zMv, zMe, zMb, zmelt, zdvo, zdva, zdM, zSs, zdMe, zdMb, zdMv47 REAL(wp) :: zMnew, zMnew1, zMnew2, zheat48 REAL(wp) :: zMbits, znMbits, zdMbitsE, zdMbitsM, zLbits, zAbits, zMbb49 REAL(wp) :: zxi, zyj, zff, z1_rday, z1_e1e2, zdt, z1_dt, z1_dt_e1e250 INTEGER :: ii, ij51 TYPE(iceberg) 52 TYPE(point) 43 INTEGER, INTENT(in) :: kt ! timestep number, just passed to print_berg 44 ! 45 INTEGER :: ii, ij 46 REAL(wp) :: zM, zT, zW, zL, zSST, zVol, zLn, zWn, zTn, znVol, zIC, zDn 47 REAL(wp) :: zMv, zMe, zMb, zmelt, zdvo, zdva, zdM, zSs, zdMe, zdMb, zdMv 48 REAL(wp) :: zMnew, zMnew1, zMnew2, zheat 49 REAL(wp) :: zMbits, znMbits, zdMbitsE, zdMbitsM, zLbits, zAbits, zMbb 50 REAL(wp) :: zxi, zyj, zff, z1_rday, z1_e1e2, zdt, z1_dt, z1_dt_e1e2 51 TYPE(iceberg), POINTER :: this, next 52 TYPE(point) , POINTER :: pt 53 53 !!---------------------------------------------------------------------- 54 54 ! … … 192 192 ! 193 193 ELSE ! Diagnose mass distribution on grid 194 z1_e1e2 194 z1_e1e2 = 1._wp / e1e2t(ii,ij) * this%mass_scaling 195 195 CALL size_budget( ii, ij, zWn, zLn, zAbits, & 196 196 & this%mass_scaling, zMnew, znMbits, z1_e1e2) … … 206 206 emp (:,:) = emp (:,:) - berg_grid%floating_melt(:,:) 207 207 emps(:,:) = emps(:,:) - berg_grid%floating_melt(:,:) 208 !! qns (:,:) = qns (:,:) + berg_grid%calving_hflx (:,:) !! heat flux not yet properly coded208 !! qns (:,:) = qns (:,:) + berg_grid%calving_hflx (:,:) !!gm heat flux not yet properly coded ==>> need it, SOLVE that! 209 209 ENDIF 210 210 ! 211 END SUBROUTINE thermodynamics 212 211 END SUBROUTINE thermodynamics 212 213 !!====================================================================== 213 214 END MODULE icbthm -
branches/2012/dev_r3337_NOCS10_ICB/NEMOGCM/NEMO/OPA_SRC/ICB/icbtrj.F90
r3359 r3370 19 19 USE in_out_manager ! NEMO IO, numout in particular 20 20 USE netcdf 21 21 ! 22 22 USE icb_oce ! define iceberg arrays 23 23 USE icbutl ! iceberg utility routines … … 26 26 PRIVATE 27 27 28 INTEGER, PRIVATE :: num_traj 29 INTEGER, PRIVATE :: n_dim, m_dim 30 INTEGER, PRIVATE :: ntrajid 31 INTEGER, PRIVATE :: numberid, nstepid, nscaling_id 32 INTEGER, PRIVATE :: nlonid, nlatid, nxid, nyid, nuvelid, nvvelid, nmassid 33 INTEGER, PRIVATE :: nuoid, nvoid, nuaid, nvaid, nuiid, nviid 34 INTEGER, PRIVATE :: nsshxid, nsshyid, nsstid, ncntid, nthkid 35 INTEGER, PRIVATE :: nthicknessid, nwidthid, nlengthid 36 INTEGER, PRIVATE :: nyearid, ndayid 37 INTEGER, PRIVATE :: nmass_of_bits_id, nheat_density_id 38 39 PUBLIC traj_init ! routine called in icbini.F90 module 40 PUBLIC traj_write ! routine called in icbrun.F90 module 41 PUBLIC traj_sync ! routine called in icbrun.F90 module 42 PUBLIC traj_end ! routine called in icbrun.F90 module 43 28 PUBLIC traj_init ! routine called in icbini.F90 module 29 PUBLIC traj_write ! routine called in icbrun.F90 module 30 PUBLIC traj_sync ! routine called in icbrun.F90 module 31 PUBLIC traj_end ! routine called in icbrun.F90 module 32 33 INTEGER :: num_traj 34 INTEGER :: n_dim, m_dim 35 INTEGER :: ntrajid 36 INTEGER :: numberid, nstepid, nscaling_id 37 INTEGER :: nlonid, nlatid, nxid, nyid, nuvelid, nvvelid, nmassid 38 INTEGER :: nuoid, nvoid, nuaid, nvaid, nuiid, nviid 39 INTEGER :: nsshxid, nsshyid, nsstid, ncntid, nthkid 40 INTEGER :: nthicknessid, nwidthid, nlengthid 41 INTEGER :: nyearid, ndayid 42 INTEGER :: nmass_of_bits_id, nheat_density_id 43 44 !!---------------------------------------------------------------------- 45 !! NEMO/OPA 3.3 , NEMO Consortium (2011) 46 !! $Id:$ 47 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 48 !!---------------------------------------------------------------------- 44 49 CONTAINS 45 50 … … 56 61 TYPE(point) , POINTER :: pt 57 62 58 IF( lk_mpp ) THEN 59 WRITE(cl_filename,'("trajectory_icebergs_",I6.6,"_",I4.4,".nc")') ktend, narea-1 60 ELSE 61 WRITE(cl_filename,'("trajectory_icebergs_",I6.6,".nc")') ktend 63 IF( lk_mpp ) THEN ; WRITE(cl_filename,'("trajectory_icebergs_",I6.6,"_",I4.4,".nc")') ktend, narea-1 64 ELSE ; WRITE(cl_filename,'("trajectory_icebergs_",I6.6 ,".nc")') ktend 62 65 ENDIF 63 66 IF ( lwp .AND. nn_verbose_level >= 0) WRITE(numout,'(2a)') 'icebergs, traj_init: creating ',TRIM(cl_filename) … … 162 165 ! End define mode 163 166 iret = NF90_ENDDEF(ntrajid) 164 167 ! 165 168 END SUBROUTINE traj_init 166 169 167 !!-------------------------------------------------------------------------168 170 169 171 SUBROUTINE traj_write( kt ) 170 !! for the moment write out each snapshot of positions 171 !! later can rewrite so that it is buffered and written out more efficiently 172 173 INTEGER, INTENT( in ) :: kt 174 ! Local variables 175 INTEGER :: iret, jn 176 CHARACTER(len=80) :: cl_filename 177 TYPE(iceberg), POINTER :: this 178 TYPE(point) , POINTER :: pt 172 !!---------------------------------------------------------------------- 173 !! *** ROUTINE traj_write *** 174 !! 175 !! ** Purpose : compute the iceberg thermodynamics. 176 !! 177 !! ** Method : - for the moment write out each snapshot of positions later 178 !! can rewrite so that it is buffered and written out more efficiently 179 !!---------------------------------------------------------------------- 180 INTEGER, INTENT( in ) :: kt 181 ! 182 INTEGER :: iret, jn 183 CHARACTER(len=80) :: cl_filename 184 TYPE(iceberg), POINTER :: this 185 TYPE(point) , POINTER :: pt 186 !!---------------------------------------------------------------------- 179 187 180 188 ! Write variables … … 218 226 219 227 this=>this%next 220 END DO228 END DO 221 229 IF( lwp .and. nn_verbose_level > 0 ) WRITE(numout,*) 'trajectory write to frame ', jn 222 230 num_traj = jn 223 231 ! 224 232 END SUBROUTINE traj_write 225 233 … … 227 235 228 236 SUBROUTINE traj_sync() 229 ! Local variables 237 !!---------------------------------------------------------------------- 238 !! *** ROUTINE traj_sync *** 239 !! 240 !! ** Purpose : 241 !!---------------------------------------------------------------------- 230 242 INTEGER :: iret 231 243 !!---------------------------------------------------------------------- 232 244 ! flush to file 233 245 iret = NF90_SYNC(ntrajid) 234 IF (iret .NE. NF90_NOERR) CALL ctl_stop('icebergs, traj_sync: nf_sync failed')235 246 IF(iret /= NF90_NOERR) CALL ctl_stop( 'icebergs, traj_sync: nf_sync failed' ) 247 ! 236 248 END SUBROUTINE traj_sync 237 249 238 !!-------------------------------------------------------------------------239 250 240 251 SUBROUTINE traj_end() 241 252 ! Local variables 242 253 INTEGER :: iret 243 254 !!---------------------------------------------------------------------- 244 255 ! Finish up 245 256 iret = NF90_CLOSE(ntrajid) 246 IF (iret .NE. NF90_NOERR) CALL ctl_stop('icebergs, traj_end: nf_close failed')247 257 IF (iret /= NF90_NOERR) CALL ctl_stop( 'icebergs, traj_end: nf_close failed' ) 258 ! 248 259 END SUBROUTINE traj_end 249 260 -
branches/2012/dev_r3337_NOCS10_ICB/NEMOGCM/NEMO/OPA_SRC/ICB/icbutl.F90
r3359 r3370 3 3 !!====================================================================== 4 4 !! *** MODULE icbutl *** 5 !! Ocean physics: various iceberg utility routines5 !! Icebergs: various iceberg utility routines 6 6 !!====================================================================== 7 7 !! History : 3.3.1 ! 2010-01 (Martin&Adcroft) Original code … … 54 54 PRIVATE create_iceberg 55 55 56 !!---------------------------------------------------------------------- 57 !! NEMO/OPA 3.3 , NEMO Consortium (2011) 58 !! $Id:$ 59 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 56 60 !!------------------------------------------------------------------------- 57 61 … … 104 108 ssh_e(jpi+1,jpj+1) = ssh_e(jpi,jpj) 105 109 CALL lbc_lnk_e( ssh_e, 'T', +1._wp, 1, 1 ) 106 110 ! 107 111 END SUBROUTINE copy_flds 108 112 109 !!-------------------------------------------------------------------------110 113 111 114 SUBROUTINE interp_flds( pi, pe1, puo, pui, pua, pssh_i, & 112 115 & pj, pe2, pvo, pvi, pva, pssh_j, & 113 & psst, pcn, phi, pff)116 & psst, pcn, phi, pff ) 114 117 !!---------------------------------------------------------------------- 115 118 !! *** ROUTINE interp_flds *** … … 138 141 pe2 = bilin_e( e2t, e2u, e2v, e2f, pi, pj ) 139 142 ! 140 puo = bilin( uo_e, pi, pj, 'U', 1, 1 )! ocean velocities141 pvo = bilin( vo_e, pi, pj, 'V', 1, 1 )143 puo = bilin( uo_e, pi, pj, 'U', 1, 1 ) ! ocean velocities 144 pvo = bilin( vo_e, pi, pj, 'V', 1, 1 ) 142 145 psst = bilin( sst_m, pi, pj, 'T', 0, 0 ) ! SST 143 146 pcn = bilin( fr_i , pi, pj, 'T', 0, 0 ) ! ice concentration 144 147 pff = bilin( ff_e , pi, pj, 'F', 1, 1 ) ! Coriolis parameter 145 148 ! 146 pua = bilin( ua_e , pi, pj, 'U', 1, 1 )! 10m wind147 pva = bilin( va_e , pi, pj, 'V', 1, 1 )! here (ua,va) are stress => rough conversion from stress to speed149 pua = bilin( ua_e , pi, pj, 'U', 1, 1 ) ! 10m wind 150 pva = bilin( va_e , pi, pj, 'V', 1, 1 ) ! here (ua,va) are stress => rough conversion from stress to speed 148 151 zcd = 1.22_wp * 1.5e-3_wp ! air density * drag coefficient 149 152 zmod = 1._wp / MAX( 1.e-20, SQRT( zcd * SQRT( pua*pua + pva*pva) ) ) 150 pua = pua * zmod! note: stress module=0 necessarly implies ua=va=0151 pva = pva * zmod153 pua = pua * zmod ! note: stress module=0 necessarly implies ua=va=0 154 pva = pva * zmod 152 155 153 156 #if defined key_lim2 || defined key_lim3 … … 167 170 END SUBROUTINE interp_flds 168 171 169 !!-------------------------------------------------------------------------170 172 171 173 REAL(wp) FUNCTION bilin( pfld, pi, pj, cd_type, kdi, kdj ) … … 223 225 END FUNCTION bilin 224 226 225 !!-------------------------------------------------------------------------226 227 227 228 REAL(wp) FUNCTION bilin_x( pfld, pi, pj ) … … 252 253 zj = pj - REAL(ij,wp) 253 254 ! 254 ! find position in this processor 255 ! find position in this processor !!gm use here mig, mjg arrays 255 256 ii = ii - nimpp + 1 256 257 ij = ij - njmpp + 1 … … 268 269 END FUNCTION bilin_x 269 270 270 !!-------------------------------------------------------------------------271 271 272 272 REAL(wp) FUNCTION bilin_e( pet, peu, pev, pef, pi, pj ) … … 299 299 ! 0.5 <= zi < 1 , 0.5 <= zj < 1 --> SW quadrant 300 300 301 zi = pi - REAL(ii,wp) 301 zi = pi - REAL(ii,wp) !!gm use here mig, mjg arrays 302 302 zj = pj - REAL(ij,wp) 303 303 … … 341 341 END FUNCTION bilin_e 342 342 343 !!-------------------------------------------------------------------------344 345 343 346 344 SUBROUTINE add_new_berg_to_list( bergvals, ptvals ) … … 365 363 END SUBROUTINE add_new_berg_to_list 366 364 367 !!------------------------------------------------------------------------- 368 369 SUBROUTINE create_iceberg(berg, bergvals, ptvals ) 370 ! Arguments 371 TYPE(iceberg), INTENT(in) :: bergvals 372 TYPE(point) , INTENT(in) :: ptvals 373 TYPE(iceberg), POINTER :: berg 374 ! local variables 375 TYPE(point) , POINTER :: pt 376 INTEGER :: istat 377 378 IF ( ASSOCIATED(berg) ) then 379 CALL ctl_stop( 'icebergs, create_iceberg: berg already associated' ) 380 ENDIF 365 366 SUBROUTINE create_iceberg( berg, bergvals, ptvals ) 367 !!---------------------------------------------------------------------- 368 !! *** ROUTINE add_new_berg_to_list *** 369 !! 370 !! ** Purpose : add a new berg to the iceberg list 371 !! 372 !! ** method : - ? 373 !!---------------------------------------------------------------------- 374 TYPE(iceberg), INTENT(in) :: bergvals 375 TYPE(point) , INTENT(in) :: ptvals 376 TYPE(iceberg), POINTER :: berg 377 ! 378 TYPE(point) , POINTER :: pt 379 INTEGER :: istat 380 !!---------------------------------------------------------------------- 381 ! 382 IF( ASSOCIATED(berg) ) CALL ctl_stop( 'icebergs, create_iceberg: berg already associated' ) 381 383 ALLOCATE(berg, STAT=istat) 382 384 IF( istat /= 0 ) CALL ctl_stop( 'failed to allocate iceberg' ) … … 385 387 berg%prev => NULL() 386 388 berg%next => NULL() 387 389 ! 388 390 ALLOCATE(pt, STAT=istat) 389 391 IF( istat /= 0 ) CALL ctl_stop( 'failed to allocate first iceberg point' ) 390 392 pt = ptvals 391 393 berg%current_point => pt 392 394 ! 393 395 END SUBROUTINE create_iceberg 394 396 395 !!-------------------------------------------------------------------------396 397 397 398 SUBROUTINE insert_berg_into_list( newberg ) … … 407 408 TYPE(iceberg), POINTER :: this, prev, last 408 409 !!---------------------------------------------------------------------- 409 410 ! 410 411 IF( ASSOCIATED( first_berg ) ) THEN 411 412 ! last = last_berg() … … 423 424 END SUBROUTINE insert_berg_into_list 424 425 425 !!-------------------------------------------------------------------------426 426 427 427 REAL(wp) FUNCTION yearday(kmon, kday, khr, kmin, ksec) 428 !!---------------------------------------------------------------------- 429 !! *** FUNCTION yearday *** 430 !! 431 !! ** Purpose : 432 !! 428 433 ! sga - improved but still only applies to 365 day year, need to do this properly 429 ! Arguments 430 INTEGER, intent(in) :: kmon, kday, khr, kmin, ksec 431 ! Local variables 434 ! 435 !!gm all these info are already known in daymod, no??? 436 !! 437 !! ** method : - ? 438 !!---------------------------------------------------------------------- 439 INTEGER, INTENT(in) :: kmon, kday, khr, kmin, ksec 440 ! 432 441 INTEGER, DIMENSION(12) :: imonths = (/ 0,31,28,31,30,31,30,31,31,30,31,30 /) 433 442 !!---------------------------------------------------------------------- 443 ! 434 444 yearday = FLOAT( SUM( imonths(1:kmon) ) ) 435 445 yearday = yearday + FLOAT(kday-1) + (FLOAT(khr) + (FLOAT(kmin) + FLOAT(ksec)/60.)/60.)/24. 436 446 ! 437 447 END FUNCTION yearday 438 448 … … 440 450 441 451 SUBROUTINE delete_iceberg_from_list( first, berg ) 442 ! Arguments 452 !!---------------------------------------------------------------------- 453 !! *** ROUTINE delete_iceberg_from_list *** 454 !! 455 !! ** Purpose : 456 !! 457 !!---------------------------------------------------------------------- 443 458 TYPE(iceberg), POINTER :: first, berg 444 459 !!---------------------------------------------------------------------- 445 460 ! Connect neighbors to each other 446 461 IF ( ASSOCIATED(berg%prev) ) THEN … … 450 465 ENDIF 451 466 IF (ASSOCIATED(berg%next)) berg%next%prev => berg%prev 452 467 ! 453 468 ! Bye-bye berg 454 469 CALL destroy_iceberg(berg) 455 470 ! 456 471 END SUBROUTINE delete_iceberg_from_list 457 472 458 !!-------------------------------------------------------------------------459 473 460 474 SUBROUTINE destroy_iceberg( berg ) 461 ! Arguments 475 !!---------------------------------------------------------------------- 476 !! *** ROUTINE destroy_iceberg *** 477 !! 478 !! ** Purpose : 479 !! 480 !!---------------------------------------------------------------------- 462 481 TYPE(iceberg), POINTER :: berg 463 482 !!---------------------------------------------------------------------- 483 ! 464 484 ! Remove any points 465 IF( ASSOCIATED( berg%current_point ) ) DEALLOCATE( berg%current_point )466 485 IF( ASSOCIATED( berg%current_point ) ) DEALLOCATE( berg%current_point ) 486 ! 467 487 ! Bye-bye berg 468 488 DEALLOCATE(berg) 469 489 ! 470 490 END SUBROUTINE destroy_iceberg 471 491 472 !!-------------------------------------------------------------------------473 492 474 493 SUBROUTINE track_berg( knum, cd_label, kt ) 475 ! Arguments 476 INTEGER, DIMENSION(nkounts) :: knum ! iceberg number 477 CHARACTER(len=*) :: cd_label 478 INTEGER :: kt ! timestep number 479 ! Local variables 494 !!---------------------------------------------------------------------- 495 !! *** ROUTINE track_berg *** 496 !! 497 !! ** Purpose : 498 !! 499 !!---------------------------------------------------------------------- 500 INTEGER, DIMENSION(nkounts) :: knum ! iceberg number 501 CHARACTER(len=*) :: cd_label ! 502 INTEGER :: kt ! timestep number 503 ! 480 504 TYPE(iceberg), POINTER :: this 481 505 LOGICAL :: match 482 506 INTEGER :: k 483 507 !!---------------------------------------------------------------------- 508 ! 484 509 this => first_berg 485 510 DO WHILE( ASSOCIATED(this) ) 486 511 match = .TRUE. 487 DO k =1,nkounts512 DO k = 1, nkounts 488 513 IF( this%number(k) /= knum(k) ) match = .FALSE. 489 514 END DO 490 IF( match ) CALL print_berg(this, kt)515 IF( match ) CALL print_berg(this, kt) 491 516 this => this%next 492 END DO493 517 END DO 518 ! 494 519 END SUBROUTINE track_berg 495 520 496 !!-------------------------------------------------------------------------497 521 498 522 SUBROUTINE print_berg( berg, kt ) 499 ! Arguments 523 !!---------------------------------------------------------------------- 524 !! *** ROUTINE print_berg *** 525 !! 526 !! ** Purpose : 527 !! 528 !!---------------------------------------------------------------------- 500 529 TYPE(iceberg), POINTER :: berg 501 530 TYPE(point) , POINTER :: pt 502 INTEGER :: kt ! timestep number 503 531 INTEGER :: kt ! timestep number 532 !!---------------------------------------------------------------------- 533 ! 504 534 pt => berg%current_point 505 535 WRITE(numicb, 9200) kt, berg%number(1), & … … 508 538 CALL flush( numicb ) 509 539 9200 FORMAT(5x,i5,2x,i10,6(2x,2f10.4)) 510 540 ! 511 541 END SUBROUTINE print_berg 512 542 513 !!-------------------------------------------------------------------------514 543 515 544 SUBROUTINE print_bergs( cd_label, kt ) 516 ! Arguments 545 !!---------------------------------------------------------------------- 546 !! *** ROUTINE print_bergs *** 547 !! 548 !! ** Purpose : 549 !! 550 !!---------------------------------------------------------------------- 517 551 CHARACTER(len=*) :: cd_label 518 552 INTEGER :: kt ! timestep number 519 ! Local variables553 ! 520 554 INTEGER :: ibergs, inbergs 521 555 TYPE(iceberg), POINTER :: this 522 556 !!---------------------------------------------------------------------- 557 ! 523 558 this => first_berg 524 559 IF( ASSOCIATED(this) ) THEN 525 560 WRITE(numicb,'(a," pe=(",i3,")")' ) cd_label, narea 526 561 WRITE(numicb,'(a8,4x,a6,12x,a5,15x,a7,19x,a3,17x,a5,17x,a5,17x,a5)' ) & 527 562 & 'timestep', 'number', 'xi,yj','lon,lat','u,v','uo,vo','ua,va','ui,vi' 528 563 ENDIF 529 564 DO WHILE( ASSOCIATED(this) ) 530 565 CALL print_berg(this, kt) 531 566 this => this%next 532 END DO567 END DO 533 568 ibergs = count_bergs() 534 569 inbergs = ibergs 535 IF( lk_mpp ) CALL mpp_sum(inbergs)536 IF ( ibergs .GT. 0 ) WRITE(numicb,'(a," there are",i5," bergs out of",i6," on PE ",i4)')&537 cd_label, ibergs, inbergs, narea538 570 IF( lk_mpp ) CALL mpp_sum(inbergs) 571 IF( ibergs > 0 ) WRITE(numicb,'(a," there are",i5," bergs out of",i6," on PE ",i4)') & 572 & cd_label, ibergs, inbergs, narea 573 ! 539 574 END SUBROUTINE print_bergs 540 575 541 !!-------------------------------------------------------------------------542 576 543 577 SUBROUTINE increment_kounter() 578 !!---------------------------------------------------------------------- 579 !! *** ROUTINE increment_kounter *** 580 !! 581 !! ** Purpose : 582 !! 544 583 ! Small routine for coping with very large integer values labelling icebergs 545 584 ! num_bergs is a array of integers … … 548 587 ! when this gets to the maximum allowed integer the second and subsequent members are 549 588 ! used to count how many times the member before cycles 550 551 ! Local variables552 INTEGER :: ii, ibig589 !!---------------------------------------------------------------------- 590 INTEGER :: ii, ibig 591 !!---------------------------------------------------------------------- 553 592 554 593 ibig = HUGE(num_bergs(1)) … … 567 606 num_bergs(1) = num_bergs(1) + jpnij 568 607 ENDIF 569 608 ! 570 609 END SUBROUTINE increment_kounter 571 610 572 !!------------------------------------------------------------------------- 573 574 INTEGER function count_bergs() 575 ! Local variables 611 612 INTEGER FUNCTION count_bergs() 613 !!---------------------------------------------------------------------- 614 !! *** FUNCTION count_bergs *** 615 !! 616 !! ** Purpose : 617 !!---------------------------------------------------------------------- 576 618 TYPE(iceberg), POINTER :: this 577 619 !!---------------------------------------------------------------------- 620 ! 578 621 count_bergs = 0 579 622 this => first_berg … … 581 624 count_bergs = count_bergs+1 582 625 this => this%next 583 END DO584 626 END DO 627 ! 585 628 END FUNCTION count_bergs 586 629 587 !!-------------------------------------------------------------------------588 630 589 631 REAL(wp) FUNCTION sum_mass( first, justbits, justbergs ) … … 624 666 END FUNCTION sum_mass 625 667 626 !!-------------------------------------------------------------------------627 668 628 669 REAL(wp) FUNCTION sum_heat( first, justbits, justbergs ) … … 663 704 END FUNCTION sum_heat 664 705 665 !!------------------------------------------------------------------------- 666 706 !!====================================================================== 667 707 END MODULE icbutl -
branches/2012/dev_r3337_NOCS10_ICB/NEMOGCM/NEMO/OPA_SRC/nemogcm.F90
r3340 r3370 163 163 #endif 164 164 165 IF( lk_diaobs ) CALL dia_obs_wri 166 IF( ln_icebergs ) CALL icb_end( nitend ) 165 IF( lk_diaobs ) CALL dia_obs_wri 166 ! 167 IF( ln_icebergs ) CALL icb_end( nitend ) 167 168 168 169 ! !------------------------!
Note: See TracChangeset
for help on using the changeset viewer.