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 cscpce.f90 24!> \brief Preparation of sending velocity variables for coupling between 25!> two instances of Code_Saturne via boundary faces. 26!> Received indformation will be transformed into boundary condition 27!> in subroutine \ref csc2cl. 28!------------------------------------------------------------------------------ 29 30!------------------------------------------------------------------------------ 31! Arguments 32!------------------------------------------------------------------------------ 33! mode name role 34!------------------------------------------------------------------------------ 35!> \param[in] nptdis 36!> \param[in] ivar variable number 37!> \param[in] locpts 38!> \param[in] vela variable value at time step beginning 39!> \param[in] coefav 40!> \param[in] coefbv 41!> \param[in] coopts 42!> \param[out] rvdis 43!______________________________________________________________________________ 44 45subroutine cscpce & 46 ( nptdis , ivar , & 47 locpts , & 48 vela , & 49 coefav , coefbv , & 50 coopts , rvdis ) 51 52!=============================================================================== 53! Module files 54!=============================================================================== 55 56use paramx 57use pointe 58use numvar 59use optcal 60use cstphy 61use cstnum 62use entsor 63use parall 64use period 65use cplsat 66use mesh 67use cs_c_bindings 68 69!=============================================================================== 70 71implicit none 72 73! Arguments 74 75integer ivar 76integer nptdis 77 78integer locpts(nptdis) 79 80double precision coopts(3,nptdis), rvdis(3,nptdis) 81double precision coefav(3 ,nfabor) 82double precision coefbv(3,3,nfabor) 83double precision vela(3,ncelet) 84 85! Local variables 86 87integer ipt , iel , isou , f_id 88integer inc , iccocg , nswrgp 89integer iwarnp , imrgrp, imligp 90 91double precision epsrgp , climgp 92double precision dx , dy , dz 93 94double precision, dimension(:,:,:), allocatable :: gradv 95 96type(var_cal_opt) :: vcopt 97 98!=============================================================================== 99 100! Allocate a temporary array 101allocate(gradv(3,3,ncelet)) 102 103call field_get_key_struct_var_cal_opt(ivarfl(ivar), vcopt) 104 105inc = 1 106iccocg = 1 107imrgrp = vcopt%imrgra 108nswrgp = vcopt%nswrgr 109imligp = vcopt%imligr 110iwarnp = vcopt%iwarni 111epsrgp = vcopt%epsrgr 112climgp = vcopt%climgr 113 114if (ivar.le.0) then 115 f_id = -1 116else 117 f_id = ivarfl(ivar) 118endif 119 120call cgdvec & 121( f_id , imrgrp , inc , nswrgp , iwarnp , imligp , & 122 epsrgp , climgp , & 123 coefav , coefbv , vela , & 124 gradv) 125 126! --- Interpolation 127 128do ipt = 1, nptdis 129 130 iel = locpts(ipt) 131 132 dx = coopts(1,ipt) - xyzcen(1,iel) 133 dy = coopts(2,ipt) - xyzcen(2,iel) 134 dz = coopts(3,ipt) - xyzcen(3,iel) 135 136 do isou = 1, 3 137 rvdis(isou,ipt) = vela(isou,iel) + gradv(1,isou,iel)*dx & 138 + gradv(2,isou,iel)*dy & 139 + gradv(3,isou,iel)*dz 140 enddo 141 142enddo 143 144! Free memory 145deallocate(gradv) 146 147!-------- 148! FORMATS 149!-------- 150!---- 151! FIN 152!---- 153 154return 155end subroutine 156