funfom Subroutine

public subroutine funfom(fc)

Objective function evaluator for VMCON author: P J Knight, CCFE, Culham Science Centre fc : output real : value of objective function at the output point This routine evaluates the value of the objective function i.e. the (normalised) figure-of-merit, at the nvar-dimensional point of interest.

Each equation for fc gives a value of the order of unity for the sake of the numerics.

Arguments

Type IntentOptional AttributesName
real(kind=dp), intent(out) :: fc

Contents

Source Code


Source Code

  subroutine funfom(fc)

    !! Objective function evaluator for VMCON
    !! author: P J Knight, CCFE, Culham Science Centre
    !! fc : output real : value of objective function at the output point
    !! This routine evaluates the value of the objective function
    !! i.e. the (normalised) figure-of-merit, at the nvar-dimensional
    !! point of interest.
    !! <P>Each equation for <CODE>fc<CODE> gives a value of the
    !! order of unity for the sake of the numerics.
    !!
    ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

		use global_variables, only: xlabel, iscan_global
		use constants, only: nout, iotty, mfile
		use constraints, only: constraint_eqns
		use cost_variables, only: concost, cfactr, cdirt, ireactor, iavail, coe
		use current_drive_variables, only: bigq, porbitlossmw, pinjmw
		use divertor_variables, only: hldiv
		use error_handling, only: idiags, fdiags, errors_on, report_error
		use heat_transport_variables, only: pnetelmw
    use numerics, only: minmax
		use physics_variables, only: fusion_power, bt, rmajor, wallmw, aspect, pohmmw
		use pf_power_variables, only: srcktpm
		use process_output, only: int_to_string3
		use tfcoil_variables, only: tfcmw
		use times_variables, only: tburn
    implicit none

    !  Arguments

    real(dp), intent(out) :: fc
!    real(c_double), intent(out) :: fc
    !  Local variables

    integer :: iab
    real(dp) :: sgn

!        write(*,*) 'Figure of merit 2 (fusion power / input power) is not used.'
    ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

    iab = abs(minmax)
    sgn = sign(1.0D0, real(minmax, kind(1.0D0)))

    !  If sgn is -1 the value of fc will be maximised
    !  If sgn is +1 the value of fc will be minimised

    select case (iab)

    case (1)  !  major radius
       fc = sgn * 0.2D0 * rmajor

    case (2)  !  fusion power / input power
        write(*,*) 'Figure of merit 2 (fusion power / input power) is not used.'
        write(*,*) 'Figure of merit 5 (fusion gain Q) is available.'
        stop 1
       ! fc = sgn * fusion_power / (pinjmw + porbitlossmw + tfcpmw + ppump/1.0D6)

    case (3)  !  neutron wall load
       fc = sgn * wallmw

    case (4)  !  TF coil + PF coil power
       fc = sgn * (tfcmw + 1.0D-3*srcktpm)/10.0D0

   case (5)  !  Q = fusion gain  Issue #540
       fc = sgn * fusion_power / (pinjmw + porbitlossmw + pohmmw)
       !fc = sgn * fusion_power / pinjmw

    case (6)  !  cost of electricity
       fc = sgn * coe/100.0D0

    case (7)  !  direct/constructed/capital cost
       if (ireactor == 0) then
          fc = sgn * cdirt/1.0D3
       else
          fc = sgn * concost/1.0D4
       end if

    case (8)  !  aspect ratio
       fc = sgn * aspect

    case (9)  !  divertor heat load
       fc = sgn * hldiv

    case (10)  !  toroidal field on axis
       fc = sgn * bt

    case (11)  !  injection power
       fc = sgn * pinjmw

    case (12)  !  hydrogen production capital cost
       ! #506 OBSOLETE
       write(*,*) 'Figure of Merit 13 (Hydrogen production) is no longer supported.'
       stop 1
    case (13)  !  hydrogen production rate
       ! #506 OBSOLETE
       write(*,*) 'Figure of Merit 13 (Hydrogen production) is no longer supported.'
       stop 1

    case (14)  !  pulse length
       fc = sgn * tburn / 2.0D4

    case (15)  !  plant availability factor (N.B. requires iavail = 1)

       if (iavail /= 1) call report_error(23)

       fc = sgn * cfactr

    case (16)  !  major radius/burn time
       fc = sgn * ( 0.95d0 * (rmajor/9.0d0) - 0.05d0 * (tburn/7200.d0) )

    case (17)  !  net electrical output
       fc = sgn * pnetelmw / 500.0d0

   case (18)  !  Null figure of merit
      fc = 1d0

   case (19)  !  major radius/burn time
      fc = sgn * ( -0.5d0 * (bigq/20.0D0) - 0.5d0 * (tburn/7200.d0) )

    case default
       idiags(1) = iab ; call report_error(24)

    end select

    !  Crude method of catching NaN errors

    if ((abs(fc) > 9.99D99).or.(fc /= fc)) then
       idiags(1) = iab ; call report_error(25)
    end if

  end subroutine funfom


© 2024

Documentation generated by FORD

PROCESS was developed by UKAEA