Jorge Gonzalez
bc8f205287
Trying to reduce warnings and unused variables in the code. This should not be in this branch.
160 lines
3.8 KiB
Fortran
160 lines
3.8 KiB
Fortran
!Module to solve the electromagnetic field
|
|
MODULE moduleEM
|
|
IMPLICIT NONE
|
|
|
|
TYPE:: boundaryEM
|
|
CHARACTER(:), ALLOCATABLE:: typeEM
|
|
INTEGER:: physicalSurface
|
|
REAL(8):: potential
|
|
|
|
CONTAINS
|
|
PROCEDURE, PASS:: apply
|
|
|
|
END TYPE boundaryEM
|
|
|
|
INTEGER:: nBoundaryEM
|
|
TYPE(boundaryEM), ALLOCATABLE:: boundEM(:)
|
|
|
|
!Information of charge and reference parameters for rho vector
|
|
REAL(8), ALLOCATABLE:: qSpecies(:)
|
|
|
|
CONTAINS
|
|
!Apply boundary conditions to the K matrix for Poisson's equation
|
|
SUBROUTINE apply(self, edge)
|
|
USE moduleMesh
|
|
IMPLICIT NONE
|
|
|
|
CLASS(boundaryEM), INTENT(in):: self
|
|
CLASS(meshEdge):: edge
|
|
INTEGER:: nNodes
|
|
INTEGER, ALLOCATABLE:: nodes(:)
|
|
INTEGER:: n
|
|
|
|
nNodes = 1
|
|
nNodes = edge%nNodes
|
|
nodes = edge%getNodes(nNodes)
|
|
|
|
DO n = 1, nNodes
|
|
SELECT CASE(self%typeEM)
|
|
CASE ("dirichlet")
|
|
mesh%K(nodes(n), :) = 0.D0
|
|
mesh%K(nodes(n), nodes(n)) = 1.D0
|
|
|
|
mesh%nodes(nodes(n))%obj%emData%type = self%typeEM
|
|
mesh%nodes(nodes(n))%obj%emData%phi = self%potential
|
|
|
|
END SELECT
|
|
|
|
END DO
|
|
|
|
END SUBROUTINE
|
|
|
|
!Assemble the source vector based on the charge density to solve Poisson's equation
|
|
SUBROUTINE assembleSourceVector(vectorF)
|
|
USE moduleMesh
|
|
USE moduleRefParam
|
|
IMPLICIT NONE
|
|
|
|
REAL(8), INTENT(out):: vectorF(1:mesh%numNodes)
|
|
REAL(8), ALLOCATABLE:: localF(:)
|
|
INTEGER, ALLOCATABLE:: nodes(:)
|
|
REAL(8), ALLOCATABLE:: rho(:)
|
|
INTEGER:: nNodes
|
|
INTEGER:: e, i, ni
|
|
CLASS(meshNode), POINTER:: node
|
|
|
|
!$OMP SINGLE
|
|
vectorF = 0.D0
|
|
!$OMP END SINGLE
|
|
|
|
!$OMP DO REDUCTION(+:vectorF)
|
|
DO e = 1, mesh%numCells
|
|
nNodes = mesh%cells(e)%obj%nNodes
|
|
nodes = mesh%cells(e)%obj%getNodes(nNodes)
|
|
!Calculates charge density (rho) in element nodes
|
|
ALLOCATE(rho(1:nNodes))
|
|
rho = 0.D0
|
|
DO i = 1, nNodes
|
|
ni = nodes(i)
|
|
node => mesh%nodes(ni)%obj
|
|
rho(i) = DOT_PRODUCT(qSpecies(:), node%output(:)%den/(vol_ref*node%v*n_ref))
|
|
|
|
END DO
|
|
|
|
!Calculates local F vector
|
|
localF = mesh%cells(e)%obj%elemF(nNodes, rho)
|
|
|
|
!Assign local F to global F
|
|
DO i = 1, nNodes
|
|
ni = nodes(i)
|
|
vectorF(ni) = vectorF(ni) + localF(i)
|
|
|
|
END DO
|
|
|
|
DEALLOCATE(localF)
|
|
DEALLOCATE(nodes, rho)
|
|
|
|
END DO
|
|
!$OMP END DO
|
|
|
|
!Apply boundary conditions
|
|
!$OMP DO
|
|
DO i = 1, mesh%numNodes
|
|
node => mesh%nodes(i)%obj
|
|
|
|
SELECT CASE(node%emData%type)
|
|
CASE ("dirichlet")
|
|
vectorF(i) = node%emData%phi
|
|
|
|
END SELECT
|
|
|
|
END DO
|
|
!$OMP END DO
|
|
|
|
END SUBROUTINE assembleSourceVector
|
|
|
|
!Solving the Poisson equation for electrostatic potential
|
|
SUBROUTINE solveElecField()
|
|
USE moduleMesh
|
|
USE moduleErrors
|
|
IMPLICIT NONE
|
|
|
|
INTEGER, SAVE:: INFO
|
|
INTEGER:: n
|
|
REAL(8), ALLOCATABLE, SAVE:: tempF(:)
|
|
EXTERNAL:: dgetrs
|
|
|
|
!$OMP SINGLE
|
|
ALLOCATE(tempF(1:mesh%numNodes))
|
|
!$OMP END SINGLE
|
|
|
|
CALL assembleSourceVector(tempF)
|
|
|
|
!$OMP SINGLE
|
|
CALL dgetrs('N', mesh%numNodes, 1, mesh%K, mesh%numNodes, &
|
|
mesh%IPIV, tempF, mesh%numNodes, info)
|
|
!$OMP END SINGLE
|
|
|
|
IF (info == 0) THEN
|
|
!Suscessful resolution of Poission equation
|
|
!$OMP DO
|
|
DO n = 1, mesh%numNodes
|
|
mesh%nodes(n)%obj%emData%phi = tempF(n)
|
|
|
|
END DO
|
|
!$OMP END DO
|
|
|
|
ELSE
|
|
!$OMP SINGLE
|
|
CALL criticalError('Poisson equation failed', 'solveElecField')
|
|
!$OMP END SINGLE
|
|
|
|
END IF
|
|
|
|
!$OMP SINGLE
|
|
DEALLOCATE(tempF)
|
|
!$OMP END SINGLE
|
|
|
|
END SUBROUTINE solveElecField
|
|
|
|
END MODULE moduleEM
|