JGonzalez/fpakc
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

Finally staying in Ampere. The other one is not quite stable. Have to figure out a way to justify this in the paper (not too hard). I have rescaled all the nodes values by the node volume to avoid issues

Browse source
This commit is contained in:
Jorge Gonzalez 2026-05-07 11:14:16 +02:00
commit 8edfd803e7
2 changed files with 13 additions and 23 deletions
Download patch file Download diff file Expand all files Collapse all files

28
src/modules/mesh/moduleMesh@boundaryEM.f90
View file

@ -239,6 +239,7 @@ submodule(moduleMesh) boundaryEM
! Update
subroutine updateFloating(self)
use moduleCaseParam, only: tauMin
use moduleRefParam, only: Vol_ref, n_ref
implicit none
class(boundaryEMGeneric), intent(inout):: self
@ -265,7 +266,7 @@ submodule(moduleMesh) boundaryEM
node => mesh%nodes(nodes(n))%obj
! Minus sign to get the flux exiting the edge
mom_nodes(n) = - dot_product(node%output(s)%mom, edge%normal)
mom_nodes(n) = - dot_product(node%output(s)%mom, edge%normal)/(node%v*Vol_ref*n_ref)
end do
@ -387,51 +388,38 @@ submodule(moduleMesh) boundaryEM
integer:: e, n, s
integer, allocatable:: nodes(:)
real(8), allocatable:: mom_nodes(:)
real(8), allocatable:: den_nodes(:)
class(meshNode), pointer:: node
real(8):: mom_center, den_center, charge_center
real(8):: currentDensity_center
select type(self)
type is(boundaryEMFreeCurrent)
self%counter = self%counter + 1
do e=1, self%nEdges
mom_center = 0.0d0
den_center = 0.0d0
charge_center = 0.0d0
currentDensity_center = 0.0d0
associate(edge => self%edges(e)%obj)
nodes = edge%getNodes(edge%nNodes)
allocate(mom_nodes(1:edge%nNodes))
allocate(den_nodes(1:edge%nNodes))
do s = 1, nSpecies
do n = 1, self%nNodes
node => mesh%nodes(nodes(n))%obj
! Minus sign as we look at the values exiting the surface (opposite to the normal)
! Density and momentum are reescaled as 'den' has no data about the node volume
! Momentum is reescaled as value at the node has no data about the node volume
mom_nodes(n) = - dot_product(node%output(s)%mom, edge%normal)/(node%v*Vol_ref*n_ref)
! den_nodes(n) = node%output(s)%den/(node%v*Vol_ref*n_ref)
end do
! mom_center = mom_center + edge%gatherF(edge%centerXi(), edge%nNodes, mom_nodes)
! den_center = den_center + edge%gatherF(edge%centerXi(), edge%nNodes, den_nodes)
! charge_center = charge_center + qSpecies(s)*edge%gatherF(edge%centerXi(), edge%nNodes, den_nodes)
mom_center = mom_center + qSpecies(s)*edge%gatherF(edge%centerXi(), edge%nNodes, mom_nodes)
currentDensity_center = currentDensity_center + qSpecies(s)*edge%gatherF(edge%centerXi(), edge%nNodes, mom_nodes)
end do
! if (den_center > 1.0e-12) then
! self%deltaElectricField(e) = self%deltaElectricField(e) + mom_center/den_center * charge_center * tauMin
self%deltaElectricField(e) = self%deltaElectricField(e) + currentDensity_center * tauMin
! end if
self%deltaElectricField(e) = self%deltaElectricField(e) + mom_center * tauMin
deallocate(nodes, mom_nodes, den_nodes)
deallocate(nodes, mom_nodes)
end associate

8
src/modules/mesh/moduleMesh@boundaryParticle.f90
View file

@ -403,6 +403,7 @@ submodule(moduleMesh) boundaryParticle
! Update
subroutine quasiNeutrality_update(self)
use moduleRefParam, only: Vol_ref, n_ref
implicit none
class(boundaryParticleGeneric), intent(inout):: self
@ -429,7 +430,7 @@ submodule(moduleMesh) boundaryParticle
do n = 1, edge%nNodes
node => mesh%nodes(nodes(n))%obj
density_nodes(n) = node%output(s)%den
density_nodes(n) = node%output(s)%den/(node%v*Vol_ref*n_ref)
end do
@ -541,6 +542,7 @@ submodule(moduleMesh) boundaryParticle
! Update
subroutine outflowAdaptive_update(self)
use moduleRefParam, only: Vol_ref, n_ref
implicit none
class(boundaryParticleGeneric), intent(inout):: self
@ -574,8 +576,8 @@ submodule(moduleMesh) boundaryParticle
do n = 1, edge%nNodes
node => mesh%nodes(nodes(n))%obj
den_nodes(n) = node%output(s)%den
mom_nodes(n) = dot_product(node%output(s)%mom, edge%normal)
den_nodes(n) = node%output(s)%den/(node%v*Vol_ref*n_ref)
mom_nodes(n) = dot_product(node%output(s)%mom, edge%normal)/(node%v*Vol_ref*n_ref)
end do