1!-------------------------------------------------------------------------------
2
3! This file is part of Code_Saturne, a general-purpose CFD tool.
4!
5! Copyright (C) 1998-2021 EDF S.A.
6!
7! This program is free software; you can redistribute it and/or modify it under
8! the terms of the GNU General Public License as published by the Free Software
9! Foundation; either version 2 of the License, or (at your option) any later
10! version.
11!
12! This program is distributed in the hope that it will be useful, but WITHOUT
13! ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
14! FOR A PARTICULAR PURPOSE.  See the GNU General Public License for more
15! details.
16!
17! You should have received a copy of the GNU General Public License along with
18! this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
19! Street, Fifth Floor, Boston, MA 02110-1301, USA.
20
21!-------------------------------------------------------------------------------
22
23!> \file clpsca.f90
24!> \brief This subroutine clips the values of a given scalar or variance.
25!>
26!-------------------------------------------------------------------------------
27
28!-------------------------------------------------------------------------------
29! Arguments
30!______________________________________________________________________________.
31!  mode           name          role                                           !
32!______________________________________________________________________________!
33!> \param[in]     iscal         scalar number
34!_______________________________________________________________________________
35
36subroutine clpsca &
37 ( iscal )
38
39!===============================================================================
40! Module files
41!===============================================================================
42
43use paramx
44use numvar
45use entsor
46use optcal
47use dimens
48use cstphy
49use cstnum
50use parall
51use field
52use cs_c_bindings
53use mesh
54
55!===============================================================================
56
57implicit none
58
59! Arguments
60
61integer          iscal
62
63! Local variables
64
65integer          ivar, iel, iflid
66integer          iclmax(1), iclmin(1), iiscav
67integer          kscmin, kscmax, f_id
68double precision vmin(1), vmax(1), vfmax
69double precision scmax, scmin
70double precision scmaxp, scminp
71
72double precision, dimension(:), pointer :: cvar_scal, cvar_scav
73
74!===============================================================================
75
76!===============================================================================
77! 1. Initialisation
78!===============================================================================
79
80! --- Numero de variable de calcul et de post associe au scalaire traite
81ivar   = isca(iscal)
82iflid  = ivarfl(ivar)
83
84! --- Numero du scalaire eventuel associe dans le cas fluctuation
85iiscav = iscavr(iscal)
86
87! Map field arrays
88call field_get_val_s(ivarfl(ivar), cvar_scal)
89
90! Key ids for clipping
91call field_get_key_id("min_scalar_clipping", kscmin)
92call field_get_key_id("max_scalar_clipping", kscmax)
93
94!===============================================================================
95! 2. Printings and clippings
96!===============================================================================
97
98! Variances are always clipped at positive values
99
100! Compute min. and max. values
101vmin(1) = cvar_scal(1)
102vmax(1) = cvar_scal(1)
103do iel = 1, ncel
104  vmin(1) = min(vmin(1),cvar_scal(iel))
105  vmax(1) = max(vmax(1),cvar_scal(iel))
106enddo
107
108! Clipping of non-variance scalars
109! And first clippings for the variances (usually 0 for min and +grand for max)
110
111iclmax(1) = 0
112iclmin(1) = 0
113
114! Get the min clipping
115call field_get_key_double(iflid, kscmin, scminp)
116call field_get_key_double(iflid, kscmax, scmaxp)
117
118if(scmaxp.gt.scminp)then
119  do iel = 1, ncel
120    if(cvar_scal(iel).gt.scmaxp)then
121      iclmax(1) = iclmax(1) + 1
122      cvar_scal(iel) = scmaxp
123    endif
124    if(cvar_scal(iel).lt.scminp)then
125      iclmin(1) = iclmin(1) + 1
126      cvar_scal(iel) = scminp
127    endif
128  enddo
129endif
130
131! Clipping of max of variances
132! Based on associated scalar values (or 0 at min.)
133if (iiscav.ne.0.and. iclvfl(iscal).eq.1) then
134
135  ! Clipping of variances
136
137
138  iclmax(1) = 0
139  iclmin(1) = 0
140
141  ! Get the corresponding scalar
142  f_id = ivarfl(isca(iiscav))
143  call field_get_val_s(f_id, cvar_scav)
144  ! Get the min clipping of the corresponding scalar
145  call field_get_key_double(f_id, kscmin, scmin)
146  call field_get_key_double(f_id, kscmax, scmax)
147
148  do iel = 1, ncel
149    vfmax = (cvar_scav(iel)-scmin)*(scmax-cvar_scav(iel))
150    if(cvar_scal(iel).gt.vfmax) then
151      iclmax(1) = iclmax(1) + 1
152      cvar_scal(iel) = vfmax
153    endif
154  enddo
155
156endif
157
158call log_iteration_clipping_field(iflid, iclmin(1), iclmax(1), vmin, vmax, iclmin(1), iclmax(1))
159
160!--------
161! Formats
162!--------
163
164!----
165! End
166!----
167
168return
169
170end subroutine
171