equilibrium_start.F90

Go to the documentation of this file.

! + + + + + + + + + + + + + + + + + + + + + + + + + + + +  
! + + + + + + + + + + + + + + + + + + + + + + + + + + + +  
MODULE equilibrium_start

CONTAINS



! + + + + + + + + + + + + + + + + + + + + + + + + + + + +  
! + + + + + + + + + + + + + + + + + + + + + + + + + + + +   

  SUBROUTINE start_profiles_consistency                               &
       (prof_flag, j0_flag, q0_flag, ext_equil,                       &
        coreprof_in, equilibrium_in, toroidfield_in, coreprof_out, equilibrium_out) 

!-------------------------------------------------------!
!     This routine checks for the consistency of        !
!     profiles (psi, q, jparallel) and equilibrium      !
!                                                       !
!     information received in: COREPROF_IN              !
!                              EQUILIBRIUM_IN.          !
!                              TOROIDFIELD_IN.          !
!                              EQUILIBRIUM_IN.          !
!                                                       !
!     information saved in:    COREPROF_OUT             !
!                              EQUILIBRIUM_OUT.         !
!                                                       !
!     controling parameter:    PROF_FLAG                !
!                              EXT_EQUIL                !
!                                                       !
!-------------------------------------------------------!
!     Source:       ---                                 !
!     Developers:   D.Kalupin                           !
!     Kontacts:     Denis.Kalupin@efda.org              !
!                                                       !
!     Comments:     created for V&V between ETS, JETTO  !
!                   and ASTRA                           !
!                                                       !
!-------------------------------------------------------!


    USE euitm_schemas
    USE euitm_routines
    USE copy_structures
    USE deallocate_structures

    USE equilibrium_input
    USE ets_wrapper_emeq_e3m
    USE ets_wrapper_bdseq
#ifdef GOT_HELENA
    USE ets_wrapper_helena
#endif
#ifdef GOT_CHEASE
  USE ets_wrapper_chease
#endif
 
    USE allocate_deallocate
    USE itm_constants
    USE equilibrium_work

    IMPLICIT NONE

    INTEGER                           :: ifail

! +++ CPO derived types:
    TYPE (type_equilibrium), POINTER  :: equilibrium_in(:)   !input CPO with geometry quantities from previous time
    TYPE (type_equilibrium), POINTER  :: equilibrium_out(:)  !output CPO with geometry quantities from previous iteration
    TYPE (type_toroidfield), POINTER  :: toroidfield_in(:)   !toroidal field, major radius, total current
    TYPE (type_coreprof),    POINTER  :: coreprof_in(:)      !input CPO with plasma profiles   
    TYPE (type_coreprof),    POINTER  :: coreprof_iter(:)    !iterration CPO with plasma profiles   
    TYPE (type_coreprof),    POINTER  :: coreprof_out(:)     !output CPO with plasma profiles   
    TYPE (type_coreprof),    POINTER  :: coreprof_new(:)     !output CPO with plasma profiles   


! +++ Dimensions:
    INTEGER                           :: prof_flag, flag     !determines primary profile
    INTEGER                           :: j0_flag             !flag for negative current, if J0_FLAG>0 negative current will be cut off
    INTEGER                           :: q0_flag             !Flag for positive dq/drho: 0-allowed, >0-cut off
    INTEGER                           :: ext_equil           !determines equilibrium module connected
    INTEGER                           :: index


! +++ Local parameters:
    REAL (R8)                         :: err_psi, err_q , err_jpar
    REAL (R8)                         :: conv

    flag  = 3                                                ! primary iterrations on current
    index = 0                                                ! Starting from parabolic profiles


!DPC
    coreprof_in(1)%rho_tor = coreprof_in(1)%rho_tor  &
         / coreprof_in(1)%rho_tor(SIZE(coreprof_in(1)%rho_tor,  dim=1))  &
         * equilibrium_in(1)%profiles_1d%rho_tor(SIZE (equilibrium_in(1)%profiles_1d%psi, dim=1))
!CPD 

    ALLOCATE                 (coreprof_iter(1))
    ALLOCATE                 (coreprof_new(1))
    CALL copy_cpo(coreprof_in(1), coreprof_iter(1))
    CALL copy_cpo(coreprof_in(1), coreprof_new(1))


! +++ Replace starting profiles with parabolas:
    CALL parabolic_prof(coreprof_in,  equilibrium_in, &
                              coreprof_out, equilibrium_out)



! +++ Provide consistent input for equilibrium:
 7  CALL equil_input(coreprof_out, toroidfield_in, equilibrium_out, equilibrium_in)



! +++ Call the selected equilibrium:
    IF     (ext_equil.EQ.1) THEN
       CALL bdseq_wrapper(equilibrium_in, equilibrium_out)
    ELSE IF(ext_equil.EQ.2) THEN 
       CALL emeq_e3m_wrapper(equilibrium_in, equilibrium_out)
    ELSE IF(ext_equil.EQ.3) THEN
#ifdef GOT_HELENA
       CALL deallocate_cpo(equilibrium_out)
       CALL helena_wrapper(equilibrium_in, equilibrium_out)
#else
       WRITE(*,*) 'No HELENA available at compile time'
       stop 'No HELENA'
#endif
    ELSE IF (ext_equil.EQ.4) THEN
#ifdef GOT_CHEASE
       CALL deallocate_cpo(equilibrium_out)
       CALL chease_wrapper(equilibrium_in, equilibrium_out)
#else
       WRITE(*,*) 'No CHEASE available at compile time'
       stop 'No CHEASE'
#endif
    END IF

!DPC
    coreprof_out(1)%rho_tor = coreprof_out(1)%rho_tor  &
         / coreprof_out(1)%rho_tor(SIZE(coreprof_out(1)%rho_tor,  dim=1))  &
         * equilibrium_out(1)%profiles_1d%rho_tor(SIZE (equilibrium_out(1)%profiles_1d%psi, dim=1))
!CPD 



! +++ Readjust profiles:
    CALL readjust_profiles            &
       (3, q0_flag, coreprof_out, equilibrium_out, coreprof_in)



! +++ Check for negative currents and total current:
    CALL correct_current_prof(1, coreprof_in, equilibrium_out, coreprof_out) 



! +++ Check for convergence of profiles :
    err_psi          = maxval(abs(1.0_r8 - abs(coreprof_iter(1)%psi%value(2:) / coreprof_out(1)%psi%value(2:))))
    err_q            = maxval(abs(1.0_r8 - abs(coreprof_iter(1)%profiles1d%q%value(2:) / coreprof_out(1)%profiles1d%q%value(2:))))
    err_jpar         = maxval(abs(1.0_r8 - abs(coreprof_iter(1)%profiles1d%jtot%value(2:) / coreprof_out(1)%profiles1d%jtot%value(2:))))

    conv             = max(err_psi, err_q , err_jpar)

    write(*,'(a,1p,4(1x,g12.5))') 'START_PROFILES_CONSISTENCY convergence (inner): ', conv, err_psi, err_q, err_jpar

! +++ Fill the iterration CPO with readjusted profiles:
    CALL copy_cpo(coreprof_out(1), coreprof_iter(1))



    IF (conv.GT.0.001_r8) goto 7



! +++ Insert experimental profiles:
    CALL readjust_profiles(prof_flag, q0_flag, coreprof_new, equilibrium_out, coreprof_out) 

 

! +++ Check for convergence of profiles :
    err_psi          = maxval(abs(1.0_r8 - abs(coreprof_out(1)%psi%value(2:) / coreprof_new(1)%psi%value(2:))))
    err_q            = maxval(abs(1.0_r8 - abs(coreprof_out(1)%profiles1d%q%value(2:) / coreprof_new(1)%profiles1d%q%value(2:))))
    err_jpar         = maxval(abs(1.0_r8 - abs(coreprof_out(1)%profiles1d%jtot%value(2:) / coreprof_new(1)%profiles1d%jtot%value(2:))))

    conv             = max(err_psi, err_q , err_jpar)

    write(*,'(a,1p,4(1x,g12.5))') 'START_PROFILES_CONSISTENCY convergence (outer): ', conv, err_psi, err_q, err_jpar
    
    CALL copy_cpo(coreprof_out(1), coreprof_new(1))

    IF (conv.GT.0.001_r8) goto 7




! +++ Deallocate CPOs:
    CALL deallocate_cpo(coreprof_iter)        
    CALL deallocate_cpo(coreprof_new)        



    RETURN


  END SUBROUTINE start_profiles_consistency

! + + + + + + + + + + + + + + + + + + + + + + + + + + + +  
! + + + + + + + + + + + + + + + + + + + + + + + + + + + +   








END MODULE equilibrium_start