My Project
programmer's documentation
Electric arcs examples

Local variables to be added

integer ifac, ii, iel
integer idim
integer izone,iesp
integer ilelt, nlelt
double precision uref2, d2s3
double precision rhomoy, dhy
double precision z1 , z2
integer ipotr, ipoti, f_id, ipotva
character(len=80) :: f_name
integer, allocatable, dimension(:) :: lstelt
double precision, dimension(:), pointer :: bfpro_rom
double precision, dimension(:,:), pointer :: cvara_potva
integer :: keyvar

Initialization and finalization

Initialization and finalization is similar to that of the base examples

Example 1

For boundary faces of color 1 assign an inlet and assign a cathode for "electric" variables.

call field_get_val_s(ibrom, bfpro_rom)
!===================
call getfbr('1', nlelt, lstelt)
!==========
do ilelt = 1, nlelt
ifac = lstelt(ilelt)
itypfb(ifac) = ientre
! Zone Number (from 1 to n)
izone = 1
! Zone localization for a given face
izfppp(ifac) = izone
rcodcl(ifac,iu,1) = 0.d0
rcodcl(ifac,iv,1) = 0.d0
rcodcl(ifac,iw,1) = 0.d0
! Turbulence
if (itytur.eq.2 .or. itytur.eq.3 &
.or. iturb.eq.50 .or. iturb.eq.60 &
.or. iturb.eq.70) then
uref2 = rcodcl(ifac,iu,1)**2 &
+rcodcl(ifac,iv,1)**2 &
+rcodcl(ifac,iw,1)**2
uref2 = max(uref2,1.d-12)
! Turbulence example computed using equations valid for a pipe.
! We will be careful to specify a hydraulic diameter adapted
! to the current inlet.
! We will also be careful if necessary to use a more precise
! formula for the dynamic viscosity use in the calculation of
! the Reynolds number (especially if it is variable, it may be
! useful to take the law from 'cs_user_physical_properties'. Here, we use by default
! the 'viscl0" value.
! Regarding the density, we have acess to its value at boundary
! faces (romb) so this value is the one used here (specifically,
! it is consistent with the processing in 'cs_user_physical_properties', in case of
! variable density)
! Hydraulic diameter
dhy = 0.075d0
rhomoy = bfpro_rom(ifac)
call turbulence_bc_inlet_hyd_diam &
( ifac, uref2, dhy, rhomoy, viscl0, rcodcl )
endif
! --- Handle Scalars
! Enthalpy in J/kg (iscalt)
! On this example we impose the value of the enthalpy
! the arbitrary value of 1.d6 corresponds to a temperature of 2200 Kelvin
! for argon at atmospheric pressure (see dp_ELE)
ii = iscalt
icodcl(ifac,isca(ii)) = 1
rcodcl(ifac,isca(ii),1) = 1.d6
! Electric potential (ipotr)
! (could corresponds also to the real part of the electrical potential if
! Joule Effect by direct conduction)
! In the Cathode example (electric arc applications),
! we impose a constant value of the electrical potential which is zero,
! assuming that the potential is equal to "ipotr + an arbitrary constant"
! (What is important for electric arc is the difference between anode and
! cathode potentials)
icodcl(ifac,ipotr) = 1
rcodcl(ifac,ipotr,1) = 0.d0
! Mass fraction of the (n-1) gas mixture components
if (ngazge .gt. 1) then
do iesp=1,ngazge-1
write(f_name,'(a13,i2.2)') 'esl_fraction_',iesp
call field_get_id(trim(f_name), f_id)
call field_get_key_int(f_id, keyvar, ii)
icodcl(ifac,ii) = 1
rcodcl(ifac,ii,1) = 0.d0
enddo
endif
! Specific model for Joule effect by direct conduction:
! Imaginary part of the potentiel (ipoti) is imposed to zero
if (ippmod(ieljou).ge. 2) then
icodcl(ifac,ipoti) = 1
rcodcl(ifac,ipoti,1) = 0.d0
endif
! Specific model for Electric arc:
! Vector Potential: Zero flux by default beacuse we don't a lot about
! vector potential (what we know, is that A is equal to zero at the infinite)
! All the boundary conditions for A are zero flux, except on some chosen faces
! where we need to impose a value in order to have a stable calculation
! (well defined problem)
! These faces are chosen where we are sure that the electrical current density
! remains very low generally far from the center of the electric arc and from
! the electrodes (see above)
if (ippmod(ielarc).ge.2) then
icodcl(ifac,ipotva ) = 3
rcodcl(ifac,ipotva ,3) = 0.d0
icodcl(ifac,ipotva + 1) = 3
rcodcl(ifac,ipotva + 1,3) = 0.d0
icodcl(ifac,ipotva + 2) = 3
rcodcl(ifac,ipotva + 2,3) = 0.d0
endif
enddo

Example 2

For boundary faces of color 5 assign an free outlet and example of electrode for Joule Effect by direct conduction.

call getfbr('5', nlelt, lstelt)
!==========
do ilelt = 1, nlelt
ifac = lstelt(ilelt)
itypfb(ifac) = isolib
! Zone Number (from 1 to n)
izone = 2
! Zone location for a given face
izfppp(ifac) = izone
! --- Handle Scalars
! Enthalpy in J/kg (By default zero flux with ISOLIB)
! Nothing to do
! Mass fraction of the (n-1) gas mixture components
! (Zero flux by defaut with ISOLIB)
! Nothing to do
! Specific model for Joule Effect by direct conduction:
! If you want to make a simulation with an imposed Power PUISIM
! (you want to get PUISIM imposed in cs_user_parameters.c and PUISIM = Amp x Volt)
! you need to impose IELCOR=1 in cs_user_parameters.c
! The boundary conditions will be scaled by COEJOU coefficient
! for example the electrical potential will be multiplied bu COEJOU
! (both real and imaginary part of the electrical potential if needed)
! COEJOU is automatically defined in order that the calculated dissipated power
! by Joule effect (both real and imaginary part if needed) is equal to PUISIM
! At the beginning of the calculation, COEJOU ie equal to 1;
! COEJOU is writing and reading in the result files.
! If you don't want to calculate with by scaling,
! you can impose directly the value.
if (ippmod(ieljou).ge. 1) then
icodcl(ifac,ipotr) = 1
if (ielcor.eq.1) then
rcodcl(ifac,ipotr,1) = 500.d0*coejou
else
rcodcl(ifac,ipotr,1) = 500.d0
endif
endif
if (ippmod(ieljou).ge. 2) then
icodcl(ifac,ipoti) = 1
if (ielcor.eq.1) then
rcodcl(ifac,ipoti,1) = sqrt(3.d0)*500.d0*coejou
else
rcodcl(ifac,ipoti,1) = sqrt(3.d0)*500.d0
endif
endif
enddo

Example 3

For boundary faces of color 2 assign a free outlet and example of anode for electric arc.

call getfbr('2', nlelt, lstelt)
!==========
do ilelt = 1, nlelt
ifac = lstelt(ilelt)
itypfb(ifac) = isolib
! Zone number (from 1 to n)
izone = 3
! Zone localization for a given face
izfppp(ifac) = izone
! --- Handle scalars
! Enthalpy in J/kg (Zero flux by default with ISOLIB)
! Nothing to do
! Real component of the electrical potential
! For electric arc model,
! ======================
! * we generally calculate the "electric variables" assuming that the total
! intensity of the electrical current is imposed (COUIMP is the value of
! the imposed total current).
! In that case, you need to impose IELCOR=1 in cs_user_parameters.c
! The "electrical variables" will be scaled by COEPOT coefficient :
! for example the electrical potential will be multiplied by COEPOT,
! Joule effect will be multipied by COEPOT * COEPOT and so on (see cs_user_electric_scaling.c)
! COEJOU is defined in cs_user_electric_scaling.c : different possibilities are described in
! cs_user_electric_scaling.c, depending on the different physics you want to simulate
! (scaling from current, from power, special model for restriking ...)
! The variable DPOT is defined: it corresponds to the electrical potential
! difference between the electrodes (Anode potential - cathode Potential).
! DPOT is calculated in cs_user_electric_scaling.c. DPOT is saved at each time step, and
! for a following calculation
! DPOT is the value of the boundary condition on anode assuming that
! the cathode potential is equel to zero.
! * It is also possible to fix the value of the potential on the anode.
! (for example, 1000 Volts).
icodcl(ifac,ipotr) = 1
if (ippmod(ielarc).ge.1 .and. ielcor .eq.1) then
rcodcl(ifac,ipotr,1) = pot_diff
else
rcodcl(ifac,ipotr,1) = 1000.d0
endif
! Mass fraction of the (n-1) gas mixture components
! zero flux by default with ISOLIB
! nothing to do
! vector Potential
! zero flux by default with ISOLIB
! nothing to do
enddo

Example 4

For boundary faces of color 3 assign a wall and example of potential vector Dirichlet condition

call getfbr('3', nlelt, lstelt)
!==========
do ilelt = 1, nlelt
ifac = lstelt(ilelt)
itypfb(ifac) = iparoi
! Zone number (from 1 to n)
izone = 4
! Zone localization for a given face
izfppp(ifac) = izone
! Wall: zero flow (zero flux for pressure)
! friction for velocities (+ turbulent variables)
! zero flux for scalars
! --- Handle scalars
! Enthalpy in J/kg (Zero flux by default)
! Nothing to do
! Real component of the electrical potential
! Zero flux by default
! Nothing to do
! Specific model for Electric arc :
! ================================
! Vector potential A (Ax, Ay, Az)
! Zero flux by default because we don't a lot about vector potential
! (what we know, is that A is equal to zero at the infinite)
! All the boundary conditions for A are zero flux, except on some chosen faces
! where we need to impose a value in order to have a stable calculation
! These faces are chosen where we are sure that the electrical current density
! remains very low generally far from the center of the electric arc and from
! the electrodes:
! On the following example, we choose to impose a "dirichlet" value for the
! 3 components of A on a small zone of the boundary located near the vertical
! free outlet of the computation domain.
! In this example, the electric arc is at the center of the computational domain,
! located on z axis (near x = 0 and y = 0).
! The x (1st) and y (the 3rd) coordinates are contained between
! -2.5 cm nd 2.5 cm:
! Ax(t, x,y,z) = Ax(t-dt, x=2.5cm, y=2.5cm, z)
! Ay(t, x,y,z) = Ay(t-dt, x=2.5cm, y=2.5cm, z)
! Az(t, x,y,z) = Az(t-dt, x=2.5cm, y=2.5cm, z)
if (ippmod(ielarc).ge.2) then
if (cdgfbo(1,ifac) .le. 2.249d-2 .or. &
cdgfbo(1,ifac) .ge. 2.249d-2 .or. &
cdgfbo(3,ifac) .le. -2.249d-2 .or. &
cdgfbo(3,ifac) .ge. 2.249d-2 ) then
call field_get_val_prev_v_by_name('vec_potential', cvara_potva)
iel = ifabor(ifac)
icodcl(ifac,ipotva) = 1
rcodcl(ifac,ipotva,1) = cvara_potva(1, iel)
icodcl(ifac,ipotva + 1) = 1
rcodcl(ifac,ipotva + 1,1) = cvara_potva(2, iel)
icodcl(ifac,ipotva + 2) = 1
rcodcl(ifac,ipotva + 2,1) = cvara_potva(3, iel)
endif
endif
enddo

Example 5

For boundary faces of color 51 assign a wall and restriking model for electric arc (anode boundaray condition).

call getfbr('51', nlelt, lstelt)
!==========
do ilelt = 1, nlelt
ifac = lstelt(ilelt)
itypfb(ifac) = iparoi
! Zone number (from 1 to n)
izone = 5
! Zone localization for a given face
izfppp(ifac) = izone
! ---- Enthalpy (J/kg) :
! imposed heat transfer coefficient
ii=iscalt
icodcl(ifac,isca(ii)) = 1
rcodcl(ifac,isca(ii),1) = 2.d4
rcodcl(ifac,isca(ii),2) = 1.d5
! Real electrical potential: anode boundary condition;
! pot_diff calculated in cs_user_electric_scaling.c
icodcl(ifac,ipotr) = 1
if (ippmod(ielarc).ge.1 .and. ielcor .eq.1) then
rcodcl(ifac,ipotr,1) = pot_diff
else
rcodcl(ifac,ipotr,1) = 100.d0
endif
! Restriking modeling:
! ===================
! example to fit depending on the case, the geometry etc...
! and also in agreement with cs_user_electric_scaling.c
if (ippmod(ielarc).ge.1 .and. ielcor .eq.1) then
if (irestrike.eq.1 .and. ntcabs.le.ntdcla+30) then
z1 = restrike_point_z - 2.d-4
if (z1.le.0.d0) z1 = 0.d0
z2 = restrike_point_z + 2.d-4
if (z2.ge.2.d-2) z2 = 2.d-2
if (cdgfbo(3,ifac).ge.z1 .and. cdgfbo(3,ifac).le.z2) then
icodcl(ifac,isca(ii)) = 1
rcodcl(ifac,isca(ii),1) = pot_diff
else
icodcl(ifac,isca(ii)) = 3
rcodcl(ifac,isca(ii),3) = 0.d0
endif
endif
endif
! Vector potential : Zero flux
if (ippmod(ielarc).ge.2) then
icodcl(ifac,ipotva ) = 3
rcodcl(ifac,ipotva ,3) = 0.d0
icodcl(ifac,ipotva + 1) = 3
rcodcl(ifac,ipotva + 1,3) = 0.d0
icodcl(ifac,ipotva + 2) = 3
rcodcl(ifac,ipotva + 2,3) = 0.d0
endif
enddo

Example 6

For boundary faces of color 4 assign a symmetry.

call getfbr('4', nlelt, lstelt)
!==========
do ilelt = 1, nlelt
ifac = lstelt(ilelt)
! Symmetries
itypfb(ifac) = isymet
! Zone number (from 1 to n)
izone = 6
! Zone localization for a given face
izfppp(ifac) = izone
! For all scalars, by default a zero flux condition is assumed
! (for potentials also)
! In Joule effect direct conduction,
! we can use an anti-symetry condition for the imaginary component of the
! electrical potential depending on the electrode configuration:
if (ippmod(ieljou).ge. 2) then
icodcl(ifac,ipoti) = 1
rcodcl(ifac,ipoti,1) = 0.d0
endif
enddo
getfbr
subroutine getfbr(fstr, facnb, faces)
Build the list of boundary faces matching a criteria string.
Definition: cs_selector_f2c.f90:111