! moduleInput: Reads JSON configuration file
MODULE moduleInput
  USE json_module
  IMPLICIT NONE

  CONTAINS
    !Main routine to read the input JSON file
    SUBROUTINE readConfig(inputFile)
      USE json_module
      USE moduleErrors
      USE moduleBoundary
      USE moduleInject
      IMPLICIT NONE

      CHARACTER(:), ALLOCATABLE, INTENT(in):: inputFile
      TYPE(json_file):: config

      !Initialize the json file variable
      CALL config%initialize()

      !Loads the config file
      CALL verboseError('Loading input file...')
      CALL config%load(filename = inputFile)

      !Reads reference parameters
      CALL readReference(config)

      !Reads case parameters
      CALL readCase(config)

      !Reads output parameters
      CALL readOutput(config)

      !Read species
      CALL readSpecies(config)

      !Read interactions between species
      CALL readInteractions(config)

      !Read boundaries
      CALL readBoundary(config)

      !Read Geometry
      CALL verboseError('Reading Geometry...')
      CALL readGeometry(config)

      !Read boundary for EM field
      CALL readEMBoundary(config)

      !Read injection of particles
      CALL verboseError('Reading Interactions between species...')
      CALL readInject(config)

      !Read parallel parameters
      CALL verboseError('Reading Parallel configuration...')
      CALL readParallel(config)

    END SUBROUTINE readConfig

    !Reads the reference parameters
    SUBROUTINE readReference(config)
      USE moduleRefParam
      USE moduleConstParam
      USE moduleErrors
      USE json_module
      IMPLICIT NONE

      TYPE(json_file), INTENT(inout):: config
      LOGICAL:: found, found_r
      CHARACTER(:), ALLOCATABLE:: object

      object = 'reference'
      !Mandatory parameters that define the case and computes derived parameters
      CALL config%get(object // '.density', n_ref, found)
      IF (.NOT. found) CALL criticalError('Reference density not found','readReference')

      CALL config%get(object // '.mass', m_ref, found)
      IF (.NOT. found) CALL criticalError('Reference mass not found','readReference')

      CALL config%get(object // '.temperature', T_ref, found)
      IF (.NOT. found) CALL criticalError('Reference temperature not found','readReference')

      CALL config%get(object // '.radius', r_ref, found_r)
      
      !Derived parameters
      v_ref     = DSQRT(kb*T_ref/m_ref)  !reference velocity
      !TODO: Make this solver dependent
      IF (found_r) THEN
        sigma_ref = PI*(r_ref+r_ref)**2    !reference cross section
        L_ref = 1.D0/(sigma_ref*n_ref)     !mean free path

      ELSE
        L_ref = DSQRT(kb*T_ref*eps_0/n_ref)/qe !Debye length
        !TODO: Obtain right sigma_ref for PIC case
        sigma_ref = PI*(4.D-10)**2        !reference cross section

      END IF
      ti_ref    = L_ref/v_ref            !reference time
      Vol_ref   = L_ref**3               !reference volume
      EF_ref    = qe*n_ref*L_ref/eps_0   !reference electric field
      Volt_ref  = EF_ref*L_ref           !reference voltage

    END SUBROUTINE readReference

    !Reads the specific case parameters
    SUBROUTINE readCase(config)
      USE moduleRefParam
      USE moduleErrors
      USE moduleCaseParam
      USE moduleSolver
      USE json_module
      IMPLICIT NONE

      TYPE(json_file), INTENT(inout):: config
      LOGICAL:: found
      CHARACTER(:), ALLOCATABLE:: object
      REAL(8):: time !simulation time in [t]
      CHARACTER(:), ALLOCATABLE:: solverType

      object = 'case'
      !Time parameters
      CALL config%get(object // '.tau', tau, found)
      IF (.NOT. found) CALL criticalError('Required parameter tau not found','readCase')

      CALL config%get(object // '.time', time, found)
      IF (.NOT. found) CALL criticalError('Required parameter time not found','readCase')

      CALL config%get(object // '.solver', solverType, found)
      IF (.NOT. found) CALL criticalError('Required parameter solver not found','readCase')

      !Allocate the type of solver
      SELECT CASE(solverType)
      CASE('2DCylNeutral')
        !Allocate the solver for neutral pusher 2D Cylindrical
        ALLOCATE(solver, source=solverCylNeutral())

      CASE('2DCylCharged')
        !Allocate the solver for charged pusher 2D Cylindrical
        ALLOCATE(solver, source=solverCylCharged())

      END SELECT

      !Convert simulation time to number of iterations
      tmax = INT(time/tau)

      !Makes tau non-dimensional
      tau = tau / ti_ref

    END SUBROUTINE readCase

    !Reads configuration for the output files
    SUBROUTINE readOutput(config)
      USE moduleErrors
      USE moduleOutput
      USE json_module
      IMPLICIT NONE

      TYPE(json_file), INTENT(inout):: config
      LOGICAL:: found
      CHARACTER(:), ALLOCATABLE:: object
      CHARACTER(8)  :: date_now=''
      CHARACTER(10) :: time_now=''

      object = 'output'
      CALL config%get(object // '.path', path, found)
      CALL config%get(object // '.trigger', triggerOutput, found)
      IF (.NOT. found) THEN
        triggerOutput = 100
        CALL warningError('Using default trigger for output file of 100 iterations')

      END IF

      !Creates output folder
      !TODO: Add option for custon name output_folder
      CALL DATE_AND_TIME(date_now, time_now)
      folder = date_now(1:4) // '-' // date_now(5:6) // '-' // date_now(7:8) // '_' &
            // time_now(1:2) // '.' // time_now(3:4) // '.' // time_now(5:6)
      CALL SYSTEM('mkdir ' // path // folder )

      CALL config%get(object // '.cpuTime', timeOutput, found)
      CALL config%get(object // '.numColl', collOutput, found)
      CALL config%get(object // '.EMField', emOutput,   found)

    END SUBROUTINE readOutput

    !Reads information about the case species
    SUBROUTINE readSpecies(config)
      USE moduleSpecies
      USE moduleErrors
      USE moduleRefParam
      USE json_module
      IMPLICIT NONE

      TYPE(json_file), INTENT(inout):: config
      CHARACTER(2):: iString
      CHARACTER(:), ALLOCATABLE:: object
      CHARACTER(:), ALLOCATABLE:: speciesType
      REAL(8):: mass, charge
      LOGICAL:: found
      INTEGER:: i

      !Gets the number of species
      CALL config%info('species', found, n_children = nSpecies)
      !Zero species means critical error
      IF (nSpecies == 0) CALL criticalError("No species found", "configRead")

      ALLOCATE(species(1:nSpecies))

      !Reads information of individual species
      DO i = 1, nSpecies
        WRITE(iString, '(I2)') i
        object = 'species(' // TRIM(iString) // ')'
        CALL config%get(object // '.type', speciesType, found)
        CALL config%get(object // '.mass', mass,        found)
        mass = mass/m_ref

        !Allocate species depending on type and assign specific parameters
        SELECT CASE(speciesType)
        CASE ("neutral")
          ALLOCATE(species(i)%obj, source=speciesNeutral())

        CASE ("charged")
          CALL config%get(object // '.charge', charge,      found)
          ALLOCATE(species(i)%obj, source=speciesCharged(q  = charge, &
                                                         qm = charge/mass))

        CASE DEFAULT
          CALL warningError("Species " // speciesType // " not supported yet")

        END SELECT
        !Assign shared parameters for all species
        CALL config%get(object // '.name',   species(i)%obj%name,   found)
        CALL config%get(object // '.weight', species(i)%obj%weight, found)
        species(i)%obj%pt = i
        species(i)%obj%m  = mass

      END DO

      !Set number of particles to 0 for init state
      !TODO: In a future, this should include the particles from init states
      nPartOld = 0

    END SUBROUTINE readSpecies

    !Reads information about interactions between species
    SUBROUTINE readInteractions(config)
      USE moduleSpecies
      USE moduleCollisions
      USE json_module
      IMPLICIT NONE

      TYPE(json_file), INTENT(inout):: config
      CHARACTER(2):: iString, kString
      CHARACTER(:), ALLOCATABLE:: object
      CHARACTER(:), ALLOCATABLE:: species_i, species_j
      CHARACTER(:), ALLOCATABLE:: crossSecFile
      CHARACTER(:), ALLOCATABLE:: crossSecFilePath
      LOGICAL:: found
      INTEGER:: nInteractions, nCollisions
      INTEGER:: i, k, ij
      INTEGER:: pt_i, pt_j

      CALL initInteractionMatrix(interactionMatrix)

      CALL config%get('interactions.folderCollisions', pathCollisions, found)

      CALL config%info('interactions.collisions', found, n_children = nInteractions)
      DO i = 1, nInteractions
        WRITE(iString, '(I2)') i
        object = 'interactions.collisions(' // TRIM(iString) // ')'
        CALL config%get(object // '.species_i', species_i, found)
        pt_i = speciesName2Index(species_i)
        CALL config%get(object // '.species_j', species_j, found)
        pt_j = speciesName2Index(species_j)
        CALL config%info(object // '.crossSections', found, n_children = nCollisions)
        ij = interactionIndex(pt_i,pt_j)
        CALL interactionMatrix(ij)%init(nCollisions)

        DO k = 1, nCollisions
          WRITE (kString, '(I2)') k
          CALL config%get(object // '.crossSections(' // TRIM(kString)// ')', crossSecFile, found)
          crossSecFilePath = pathCollisions // crossSecFile
          CALL interactionMatrix(ij)%collisions(k)%obj%init(crossSecFilePath, species(pt_i)%obj%m, species(pt_j)%obj%m)

        END DO

      END DO

    END SUBROUTINE readInteractions

    !Reads boundary conditions for the mesh
    SUBROUTINE readBoundary(config)
      USE moduleBoundary
      USE moduleErrors
      USE json_module
      IMPLICIT NONE

      TYPE(json_file), INTENT(inout):: config
      CHARACTER(2):: istring
      CHARACTER(:), ALLOCATABLE:: object
      LOGICAL:: found
      INTEGER:: i

      CALL config%info('boundary', found, n_children = nBoundary)
      ALLOCATE(boundary(1:nBoundary))
      DO i = 1, nBoundary
        WRITE(istring, '(i2)') i
        object = 'boundary(' // trim(istring) // ')'

        ALLOCATE(boundaryGeneric:: boundary(i)%obj)

        CALL config%get(object // '.type',            boundary(i)%obj%boundaryType,    found)
        CALL config%get(object // '.physicalSurface', boundary(i)%obj%physicalSurface, found)
        boundary(i)%obj%id = i

      END DO

    END SUBROUTINE readBoundary

    !Read the geometry (mesh) for the case
    SUBROUTINE readGeometry(config)
      USE moduleMesh
      USE moduleMeshCylRead
      USE moduleErrors
      USE moduleOutput
      USE json_module
      IMPLICIT NONE

      TYPE(json_file), INTENT(inout):: config
      LOGICAL:: found
      CHARACTER(:), ALLOCATABLE:: geometryType, meshType, meshFile
      CHARACTER(:), ALLOCATABLE:: fullPath

      !Selects the type of geometry.
      CALL config%get('geometry.type', geometryType, found)
      SELECT CASE(geometryType)
      CASE ("2DCyl")
        !Creates a 2D cylindrical mesh
        ALLOCATE(meshCylGeneric:: mesh)

        !Gets the type of mesh
        CALL config%get('geometry.meshType', meshType, found)
        SELECT CASE(meshType)
        CASE ("gmsh")
          !Gets the gmsh file
          CALL config%get('geometry.meshFile', meshFile, found)

        CASE DEFAULT
          CALL criticalError("Mesh type " // meshType // " not supported.", "readGeometry")
          
        END SELECT
        !Reads the mesh
        fullpath = path // meshFile
        CALL mesh%readMesh(fullPath)

      CASE DEFAULT
        CALL criticalError("Geometry type " // geometryType // " not supported.", "readGeometry")

      END SELECT

    END SUBROUTINE readGeometry

    SUBROUTINE readEMBoundary(config)
      USE moduleMesh
      USE moduleOutput
      USE moduleErrors
      USE moduleEM
      USE moduleRefParam
      USE moduleSpecies
      USE json_module
      IMPLICIT NONE

      TYPE(json_file), INTENT(inout):: config
      CHARACTER(:), ALLOCATABLE:: object
      LOGICAL:: found
      CHARACTER(2):: istring
      INTEGER:: i, e, s
      INTEGER:: info
      EXTERNAL:: dgetrf

      CALL config%info('boundaryEM', found, n_children = nBoundaryEM)
      
      IF (found) ALLOCATE(boundEM(1:nBoundaryEM))
      
      DO i = 1, nBoundaryEM
        WRITE(istring, '(i2)') i
        object = 'boundaryEM(' // trim(istring) // ')'

        CALL config%get(object // '.type', boundEM(i)%typeEM, found)

        SELECT CASE(boundEM(i)%typeEM)
        CASE ("dirichlet")
          CALL config%get(object // '.potential', boundEM(i)%potential, found)
          IF (.NOT. found) &
            CALL criticalError('Required parameter "potential" for Dirichlet boundary condition not found', 'readEMBoundary')
          boundEM(i)%potential = boundEM(i)%potential/Volt_ref

          CALL config%get(object // '.physicalSurface', boundEM(i)%physicalSurface, found)
          IF (.NOT. found) &
            CALL criticalError('Required parameter "physicalSurface" for Dirichlet boundary condition not found', 'readEMBoundary')

        CASE DEFAULT
          CALL criticalError('Boundary type ' // boundEM(i)%typeEM // ' not yet supported', 'readEMBoundary')

        END SELECT

      END DO

      ALLOCATE(qSpecies(1:nSpecies))
      DO s = 1, nSpecies
        SELECT TYPE(sp => species(s)%obj)
        TYPE IS (speciesCharged)
          qSpecies(s) = sp%q

        CLASS DEFAULT
          qSpecies(s) = 0.D0

        END SELECT

      END DO

      !TODO: Improve this
      IF (ALLOCATED(boundEM)) THEN
        DO e = 1, mesh%numEdges
          IF (ANY(mesh%edges(e)%obj%physicalSurface == boundEM(:)%physicalSurface)) THEN
            DO i = 1, nBoundaryEM
              IF (mesh%edges(e)%obj%physicalSurface == boundEM(i)%physicalSurface) THEN
                CALL boundEM(i)%apply(mesh%edges(e)%obj)

              END IF
            END DO
          END IF
        END DO
      END IF

      !Compute the PLU factorization of K once boundary conditions have been read
      CALL dgetrf(mesh%numNodes, mesh%numNodes, mesh%K, mesh%numNodes, mesh%IPIV, info)
      IF (info /= 0) CALL criticalError('Factorization of K matrix failed', 'readEMBoundary')

    END SUBROUTINE readEMBoundary

    !Reads the injection of particles from the boundaries
    SUBROUTINE readInject(config)
      USE moduleSpecies
      USE moduleErrors
      USE moduleInject
      USE json_module
      IMPLICIT NONE

      TYPE(json_file), INTENT(inout):: config
      INTEGER:: i
      CHARACTER(2):: istring
      CHARACTER(:), ALLOCATABLE:: object
      LOGICAL:: found
      CHARACTER(:), ALLOCATABLE:: speciesName
      CHARACTER(:), ALLOCATABLE:: name
      REAL(8):: v
      REAL(8), ALLOCATABLE:: T(:), normal(:)
      REAL(8):: flow
      CHARACTER(:), ALLOCATABLE:: units
      INTEGER:: physicalSurface
      INTEGER:: pt

      CALL config%info('inject', found, n_children = nInject)
      ALLOCATE(inject(1:nInject))
      nPartInj = 0
      DO i = 1, nInject
        WRITE(istring, '(i2)') i
        object = 'inject(' // trim(istring) // ')'

        !Find species
        CALL config%get(object // '.species', speciesName, found)
        pt = speciesName2Index(speciesName)

        CALL config%get(object // '.name',            name,            found)
        CALL config%get(object // '.v',               v,               found)
        CALL config%get(object // '.T',               T,               found)
        CALL config%get(object // '.n',               normal,          found)
        CALL config%get(object // '.flow',            flow,            found)
        CALL config%get(object // '.units',           units,           found)
        CALL config%get(object // '.physicalSurface', physicalSurface, found)

        CALL inject(i)%init(i, v, normal, T, flow, units, pt, physicalSurface)

      END DO

      !Allocate array for injected particles
      IF (nPartInj > 0) THEN
        ALLOCATE(partInj(1:nPartInj))

      END IF
      
    END SUBROUTINE readInject

    SUBROUTINE readParallel(config)
      USE moduleParallel
      USE moduleErrors
      USE json_module
      IMPLICIT NONE

      TYPE(json_file), INTENT(inout):: config
      CHARACTER(:), ALLOCATABLE:: object
      LOGICAL:: found

      !Reads OpenMP parameters
      object = 'parallel.OpenMP'

      CALL config%get(object // '.nThreads', openMP%nThreads, found)

      IF (.NOT. found) THEN
        openMP%nthreads = 8
        CALL warningError('No OpenMP configuration detected, using 8 threads as default.')

      END IF

      CALL initParallel()

    END SUBROUTINE readParallel


END MODULE moduleInput