Subversion Repositories lagranto.arpege

module evaluate

! The user can assign values to parameters that can be used in expressions with 

! the subroutine defparam. The calling syntax is

!

!     call defparam(symbol,value) or call defparam(symbol,expr)

!

! where symbol is the desired parameter name; value is a real, integer, or 

! complex variable (single or double precision); and expr is a string

! containing an expression to be evaluated. The value obtained by evaluating the

! expression expr is associated with the parameter symbol. Parameter names must

! begin with a letter (a-z, A-Z) and must not be longer than 24 characters.

! Parameter names are not case dependent.

!

! An expression can be evaluated with the subroutine evalexpr. The calling

! syntax is 

!

!          call evalexpr(expr,value)

!

! where expr is a string containing the expression to be evaluated; value is the

! result (single or double precision real, complex or integer). The 

! expression can contain the arithmetic operations +, -, *, /, or ^ as well as

! the functions sin, cos, tan, log, ln, abs, exp, sqrt, real, imag, conjg, and 

! ang (the function ang calculates the phase angle of its complex argument). The

! expression can also contain numerical values and previously defined parameters

! Grouping by nested levels of parentheses is also allowed. The parameters pi

! and i (imaginary unit) are predefined. Complex numbers can be entered as a+i*b 

! if the parameter i has not been redefined by the user. Complex numbers can also

! be entered using complex(a,b).

! Example expression: 

!

!          conjg(((cos(x) + sqrt(a+i*b))^2+complex(ln(1.6e-4),20))/2)

!

! An equation of the form <symbol> = <expression> can be evaluated using the 

! subroutine evaleqn. The calling syntax is

!

!          call evaleqn(eqn)       

!

! where eqn is a string containing the equation. The right-hand-side of the

! equation is evaluated and assigned to the symbol given by the left-hand-side.

!

! The value assigned to a symbol can be retrieved using the subroutine getparam.

! The calling syntax is

!

!          call getparam(sym,value)

!

! where sym is a symbol string; value is a numeric variable (any of the six 

! standard types).

!

! The symbols and their values in the symbol table can be listed using the

! subroutine listvar. The variable ierr is always available following a call

! to any of the above subroutines and is zero if there were no errors. The

! possible nonzero values for ierr are

!

! 1       Expression empty

! 2       Parentheses don't match

! 3       Number string does not correspond to a valid number

! 4       Undefined symbol

! 5       Less than two operands for binary operation

! 6       No operand for unary plus or minus operators

! 7       No argument(s) for function

! 8       Zero or negative real argument for logarithm

! 9       Negative real argument for square root

! 10      Division by zero

! 11      Improper symbol format

! 12      Missing operator

! 13      Undefined function

! 14      Argument of tangent function a multiple of pi/2

!

use precision

use strings

save

private

public :: valuep,evalexpr,defparam,evaleqn,getparam,listvar,ierr

type item                         

  character(len=24):: char

  character :: type

end type item

type param                        

  character (len=24):: symbol

  complex(kc8):: value

end type param

interface defparam                

  module procedure strdef       ! value given by expression      

  module procedure valdef_dc    ! Double precision complex value

  module procedure valdef_sc    ! Single precision complex value

  module procedure valdef_dr    ! Double precision real value

  module procedure valdef_sr    ! Single precision real value

  module procedure valdef_di    ! Double precision integer value

  module procedure valdef_si    ! Single precision integer value

end interface

interface evalexpr

  module procedure evalexpr_dc  ! Double precision complex result

  module procedure evalexpr_sc  ! Single precision complex result

  module procedure evalexpr_dr  ! Double precision real result

  module procedure evalexpr_sr  ! Single precision real result

  module procedure evalexpr_di  ! Double precision integer result

  module procedure evalexpr_si  ! Single precision integer result

end interface

interface getparam

  module procedure getparam_dc  ! Double precision complex result

  module procedure getparam_sc  ! Single precision complex result

  module procedure getparam_dr  ! Double precision real result

  module procedure getparam_sr  ! Single precision real result

  module procedure getparam_di  ! Double precision integer result

  module procedure getparam_si  ! Single precision integer result

end interface

integer,parameter :: numtok=100  ! Maximum number of tokens 

type(param) :: params(100)       ! Symbol table

integer :: nparams=0,itop,ibin

complex(kc8) :: valstack(numtok) ! Stack used in evaluation of expression

type(item):: opstack(numtok)     ! Operator stack used in conversion to postfix

integer :: ierr                  ! Error flag

!**********************************************************************

contains

!**********************************************************************

SUBROUTINE EVALEXPR_DC(expr,val)    ! Evaluate expression expr for

                                    ! val double precision complex

character (len=*),intent(in) :: expr

complex(kc8) :: val

character (len=len(expr)+1) :: tempstr

character :: cop

integer :: isp(numtok)          ! On stack priority of operators in opstack

integer :: lstr

complex(kc8) :: cval,oper1,oper2

real(kr8) :: valr,vali

type(item):: token(numtok)      ! List of tokens ( a token is an operator or 

                                ! operand) in postfix order

type(item) :: x,junk,tok

ierr=0

token(1:)%char=' '

if(nparams == 0) then                  ! Initialize symbol table

  params(1)%symbol='PI'

  params(1)%value=(3.14159265358979_kr8,0.0_kr8)

  params(2)%symbol='I'

  params(2)%value=(0.0_kr8,1.0_kr8)

  nparams=2

end if

if(len_trim(expr) == 0) then           ! Expression empty

  ierr=1

  write(*,*) 'Error: expression being evaluated is empty'

  return

end if

tempstr=adjustl(expr)

call removesp(tempstr)   ! Removes spaces, tabs, and control characters

! ****************************************************************************

! STEP 1:  Convert string to token array. Each token is either an operator or

!          an operand. Token array will be in postfix (reverse Polish) order.

!*****************************************************************************

ntok=0

ibin=0

itop=0

do

  lstr=len_trim(tempstr)

  call get_next_token(tempstr(1:lstr),tok,icp,insp)

  select case(tok%type)

  case('S')

    ntok=ntok+1

    token(ntok)=tok

  case('E')

    do 

      if(itop < 1)exit

      call popop(x)        ! Output remaining operators on stack

      ntok=ntok+1

      token(ntok)=x

    end do

    ntok=ntok+1

    token(ntok)=tok

    exit

  case('R')  ! Token is right parenenthesis

    do 

      if(opstack(itop)%type == 'L') exit  ! Output operators on stack down

      call popop(x)                       ! to left parenthesis

      ntok=ntok+1

      token(ntok)=x

    end do                             

    call popop(junk)                      ! Remove left parenthesis from stack

    if(opstack(itop)%type == 'F') then    ! Output function name if present

      call popop(x)

      ntok=ntok+1

      token(ntok)=x

    end if                               

  case('D')  ! Token is comma

    do 

      if(opstack(itop)%type == 'L') exit  ! Output operators on stack down

      call popop(x)                       ! to left parenthesis

      ntok=ntok+1

      token(ntok)=x

    end do                              

  case('U','B','L','F') ! Token is operator, left parenthesis or function name

    do 

      if(isp(itop) < icp) exit            ! Output operators on stack having

      call popop(x)                       ! an instack priority that is

      ntok=ntok+1                         ! greater than or equal to the

      token(ntok)=x                       ! priority of the incoming operator  

    end do                                

    call pushop(tok)     ! Put incoming operator on stack                       

    isp(itop)=insp

  end select

end do

isum=0                                 ! Error check for matching parentheses

do i=1,ntok

  if(token(i)%type == 'L' ) isum=isum+1

  if(token(i)%type == 'R' ) isum=isum-1

end do

if(isum /= 0) then

  ierr=2

  write(*,*) 'Error in the evaluation of the expression ',trim(expr)

  write(*,*) "Parentheses don't match"

  write(*,*)

  return

end if

!*****************************************************************************

! STEP 2: Evaluate token string in postfix order

!*****************************************************************************

itop=0

do i=1,ntok

  x=token(i)

  select case(x%type)

  case('E')  ! Token is end token

    if(itop>1) then                

      ierr=12

      write(*,*) 'Error: missing operator in expression ',trim(expr)

      write(*,*)

      return

    end if

    call popval(val)               ! Final result left on stack of values

    exit

  case('S')  ! Token is operand

    call valuep(x%char,cval)       ! Evaluate operand

    if(ierr/=0) return

    call pushval(cval)             ! Put value of operand on stack

  case('B')  ! Token is a binary operator

    if(itop < 2) then

      ierr=5

      write(*,*) 'Error in evaluation of expression ',trim(expr)

      write(*,*) 'Less than two operands for binary operator  '&

                 ,trim(x%char)

      write(*,*)

      return

    end if                         

    call popval(oper1)             ! Pull off top two values from stack

    call popval(oper2)

    select case(trim(x%char))      ! Perform operation on values

    case('^')

      cval=oper2**oper1

    case('*')

      cval=oper2*oper1

    case('/')

      if(oper1 == (0._kr8,0._kr8)) then

        ierr=10

        write(*,*) 'Error in expression ',trim(expr)

        write(*,*) 'Division by zero'

        write(*,*)

        return

      end if

      cval=oper2/oper1

    case('+')

      cval=oper2+oper1

    case('-')

      cval=oper2-oper1

    end select

    call pushval(cval)             ! Put result back on stack

  case('U')  ! Token is unary operator

    if(itop == 0) then

      ierr=6

      write(*,*) 'Error in expression ',trim(expr)

      write(*,*) 'No operand for unary operator ',trim(x%char)

      write(*,*)

      return

    else

      call popval(oper1)           ! Pull top value off stack

    end if

    select case(trim(x%char))      ! Operate on value

    case('+')

      cval=oper1

    case('-')

      cval=-oper1

    end select

    call pushval(cval)             ! Put result back on stack

  case('F')  ! Token is a function name

    if(itop == 0) then

      ierr=7

      write(*,*) 'Error in expression ',trim(expr)

      write(*,*) 'Missing argument(s) for function ',trim(x%char)

      write(*,*)

      return

    else  

      call popval(oper1)           ! Pull top value off stack

    end if 

    tempstr=uppercase(x%char)

    select case(trim(tempstr))      ! Evaluate function

    case('SIN')

      cval=sin(oper1)

    case('COS')

      cval=cos(oper1)

    case('TAN')

      oper2=cos(oper1)

      if(abs(oper2) == 0.0_kr8) then

        ierr=14

        write(*,*) 'Error: argument of tan function a multiple',&

        ' of pi/2 in expression ',trim(expr)

        write(*,*)

        return

      else 

        cval=sin(oper1)/oper2

      endif

    case('SQRT')

      if(real(oper1,kr8) < 0. .and. aimag(oper1)==0.) then

        ierr=9

        write(*,*) 'Warning: square root of negative real number',&

                   ' in expression ',trim(expr)

        write(*,*)

      end if

      cval=sqrt(oper1)

    case('ABS')

      cval=abs(oper1)

    case('LN')

      if(real(oper1,kr8) <= 0. .and. aimag(oper1)==0.) then

        ierr=8

        write(*,*) 'Error: negative real or zero argument for',&

                   ' natural logarithm in expression ',trim(expr)

        write(*,*)

        return

      end if

      cval=log(oper1)

    case('LOG')

      if(real(oper1,kr8) <= 0. .and. aimag(oper1)==0.) then

        ierr=8

        write(*,*) 'Error: negative real or zero argument for base',&

                   '10 logarithm in expression ',trim(expr)

        write(*,*)

        return

      end if

      cval=log(oper1)/2.30258509299405_kr8

    case('EXP')

      cval=exp(oper1)

    case('COMPLEX')

      if(itop == 0) then

        ierr=7

        write(*,*) 'Error in expression ',trim(expr)

        write(*,*) 'Missing argument(s) for function ',trim(x%char)

        write(*,*)

        return

      else  

        call popval(oper2)  ! Pull second argument off stack

      end if 

      valr=real(oper2,kr8)

      vali=real(oper1,kr8)

      cval=cmplx(valr,vali,kc8)

    case('CONJG')

      cval=conjg(oper1)

    case('ANG')

      cval=atan2(aimag(oper1),real(oper1,kr8))

    case('REAL')

      cval=real(oper1,kr8)

    case('IMAG')

      cval=aimag(oper1)

    case default ! Undefined function

      ierr=13

      write(*,*) 'Error: the function ',trim(x%char), ' is undefined',&

                 ' in the expression ',trim(expr)

      write(*,*)

      return

    end select

    call pushval(cval)    ! Put result back on stack

  end select

end do

end subroutine evalexpr_dc

!**********************************************************************

SUBROUTINE GET_NEXT_TOKEN(str,tok,icp,isp)

character(len=*) :: str

character :: cop,chtemp

type(item) :: tok

integer :: icp

lstr=len_trim(str)

if(lstr == 0) then

  tok%char='#'             ! Output end token

  tok%type='E'

  return

end if

ipos=scan(str,'+-*/^(),')  ! Look for an arithmetic operator 

                           ! + - * / ^ ( ) or ,

cop=str(ipos:ipos)                 

select case (ipos)              

case(0)    ! Operators not present

  ntok=ntok+1

  tok%char=str

  tok%type='S'

  str=''

  icp=0

  isp=0

case(1) 

  tok%char=cop

  select case(cop)

  case('+','-')

    if(ibin==0) then

      tok%type='U'

      icp=4

      isp=3

    else

      tok%type='B'

      icp=1

      isp=1

    end if

    ibin=0

  case('*','/')

    tok%type='B'

    icp=2

    isp=2

    ibin=0

  case('^')

    tok%type='B'

    icp=4

    isp=3

    ibin=0

  case('(')

    tok%type='L'

    icp=4

    isp=0

    ibin=0

  case(')')

    tok%type='R'

    icp=0

    isp=0

    ibin=1

  case(',')

    tok%type='D'

    icp=0

    isp=0

    ibin=0

  end select

  str=str(2:)

case(2:)

  select case(cop)

  case('(')

    tok%char=str(1:ipos-1)

    tok%type='F'

    icp=4

    isp=0

    ibin=0

    str=str(ipos:)

  case('+','-')

    chtemp=uppercase(str(ipos-1:ipos-1))

    if(is_letter(str(1:1)).eqv..true. .or. chtemp/='E') then

      tok%char=str(1:ipos-1)

      tok%type='S'

      icp=0

      isp=0

      ibin=1

      str=str(ipos:)

    else

      inext=scan(str(ipos+1:),'+-*/^(),')

      if(inext==0) then

        tok%char=str

        tok%type='S'

        icp=0

        isp=0

        ibin=0

        str=''

      else

        tok%char=str(1:ipos+inext-1)

        tok%type='S'

        icp=0

        isp=0

        ibin=1

        str=str(ipos+inext:)

      end if

    end if

  case default

    tok%char=str(1:ipos-1)

    tok%type='S'

    icp=0

    isp=0

    ibin=1

    str=str(ipos:)

  end select

end select

end subroutine get_next_token

!**********************************************************************

SUBROUTINE EVALEXPR_SC(expr,val)    ! Evaluate expression expr for

                                    ! val single precision complex

character(len=*) :: expr

complex(kc4) :: val

complex(kc8) :: vald

call evalexpr_dc(expr,vald)

val=vald

end subroutine evalexpr_sc

!**********************************************************************

SUBROUTINE EVALEXPR_SR(expr,val)    ! Evaluate expression expr for

                                    ! val single precision real

character(len=*) :: expr

real(kr4) :: val

complex(kc8) :: vald

call evalexpr_dc(expr,vald)

val=real(vald)

end subroutine evalexpr_sr

!**********************************************************************

SUBROUTINE EVALEXPR_DR(expr,val)    ! Evaluate expression expr for 

                                    ! val double precision real

character(len=*) :: expr

real(kr8) :: val

complex(kc8) :: vald

call evalexpr_dc(expr,vald)

val=real(vald,kr8)

end subroutine evalexpr_dr

!**********************************************************************

SUBROUTINE EVALEXPR_SI(expr,ival)   ! Evaluate expression expr for 

                                    ! ival single precision integer

character(len=*) :: expr

integer(ki4) :: ival

complex(kc8) :: vald

call evalexpr_dc(expr,vald)

ival=nint(real(vald,kr8),ki4)

end subroutine evalexpr_si

!**********************************************************************

SUBROUTINE EVALEXPR_DI(expr,ival)   ! Evaluate expression expr for 

                                    ! ival double precision integer

character(len=*) :: expr

integer(ki8) :: ival

complex(kc8) :: vald

call evalexpr_dc(expr,vald)

ival=nint(real(vald,kr8),ki8)

end subroutine evalexpr_di

!**********************************************************************

SUBROUTINE VALDEF_DC(sym,val)    ! Associates sym with val in symbol table,

                                 ! val double precision complex

character(len=*) :: sym

character(len=len_trim(sym)) :: usym

complex(kc8) :: val

ierr=0

if(nparams == 0) then               ! Initialize symbol table

  params(1)%symbol='PI'

  params(1)%value=(3.14159265358979_kr8,0.0_kr8)

  params(2)%symbol='I'

  params(2)%value=(0.0_kr8,1.0_kr8)

  nparams=2

end if

! Assign val to sym if sym is already in symbol table

usym=uppercase(sym)

if(is_letter(sym(1:1)).eqv..false. .or. len_trim(sym)>24) then

  ierr=11

  write(*,*) 'Error: symbol ',trim(sym),' has improper format'

  write(*,*)

  return

end if

do i=1,nparams

  if(trim(usym)==trim(params(i)%symbol)) then

    params(i)%value=val

    return

  end if

end do

nparams=nparams+1    ! Otherwise assign val to new symbol sym

params(nparams)%symbol=usym

params(nparams)%value=val

end subroutine valdef_dc

!**********************************************************************

SUBROUTINE VALDEF_SC(sym,val)     ! Associates sym with val in symbol table,

                                  ! val single precision complex

character(len=*) :: sym

complex(kc4) :: val

complex(kc8) :: vald

vald=val

call valdef_dc(sym,vald)

end subroutine valdef_sc

!**********************************************************************

SUBROUTINE VALDEF_DR(sym,val)    ! Associates sym with val in symbol table,

                                 ! val double precision real

character(len=*) :: sym

real(kr8) :: val

complex(kc8) :: vald

vald=cmplx(val,0.0_kr8,kc8)

call valdef_dc(sym,vald)

end subroutine valdef_dr

!**********************************************************************

SUBROUTINE VALDEF_SR(sym,val)    ! Associates sym with val in symbol table,

                                 ! val single precision real

character(len=*) :: sym

real(kr4) :: val

complex(kc8) :: vald

vald=cmplx(val,0.0,kc8)

call valdef_dc(sym,vald)

end subroutine valdef_sr

!**********************************************************************

SUBROUTINE VALDEF_DI(sym,ival)   ! Associates sym with ival in symbol table,

                                 ! ival double precision integer 

character(len=*) :: sym

integer(ki8) :: ival

complex(kc8) :: vald

vald=cmplx(real(ival,kr8),0.0_kr8,kc8)

call valdef_dc(sym,vald)

end subroutine valdef_di

!**********************************************************************

SUBROUTINE VALDEF_SI(sym,ival)   ! Associates sym with ival in symbol table,

                                 ! ival single precision integer

character(len=*) :: sym

integer(ki4) :: ival

complex(kc8) :: vald

vald=cmplx(real(ival,kr8),0.0,kc8)

call valdef_dc(sym,vald)

end subroutine valdef_si

!**********************************************************************

SUBROUTINE STRDEF(sym,expr)      ! Associates sym with the value of the

                                 ! expression expr

character(len=*) :: sym,expr

complex(kc8) :: val

if(nparams == 0) then            ! Initialize symbol table

  params(1)%symbol='PI'

  params(1)%value=(3.14159265358979_kr8,0.0_kr8)

  params(2)%symbol='I'

  params(2)%value=(0.0_kr8,1.0_kr8)

  nparams=2

end if

call evalexpr_dc(expr,val)       ! val is value of expression expr

if(ierr==0 .or. ierr==9) then

  call valdef_dc(sym,val)          ! Assign val to symbol sym

end if

end subroutine strdef

!**********************************************************************

SUBROUTINE VALUEP(xinchar,cval)  ! Finds double precision complex value 

                                 ! corresponding to number string xinchar 

                                 ! or value in symbol table corresponding 

                                 ! to symbol name xinchar.

character (len=*):: xinchar

complex(kc8) :: cval

real(kr8) :: rval

ierr=0

if(is_letter(xinchar(1:1)).eqv..true.) then   ! xinchar is a symbol

  call getparam(xinchar,cval)

else                               ! xinchar is a number string

  call value(xinchar,rval,ios)     ! rval is the value of xinchar

  if(ios > 0) then

    ierr=3

    write(*,*) 'Error: number string ',trim(xinchar),' does not correspond to a valid number' 

    write(*,*)

  end if

  cval=cmplx(rval,0.0_kr8,kc8)

  return

end if

end subroutine valuep

!**********************************************************************

SUBROUTINE PUSHOP(op)  ! Puts an operator on operator stack

type(item):: op

itop=itop+1

if(itop > numtok) then

  write(*,*) 'Error: operator stack overflow in evaluation of expression'

  write(*,*)

  return

end if

opstack(itop)=op

end subroutine pushop

SUBROUTINE POPOP(op) ! Takes top operator of operator stack and assigns it to op

type(item):: op

op=opstack(itop)

itop=itop-1

end subroutine popop

SUBROUTINE PUSHVAL(val) ! Puts value on value stack

complex(kc8) :: val

itop=itop+1

if(itop > numtok) then

  write(*,*) 'Error: value stack overflow in evaluation of expression'

  write(*,*)

  return

end if

valstack(itop)=val

end subroutine pushval

SUBROUTINE POPVAL(val) ! Takes top value off value stack and assigns it to val

complex(kc8) :: val

val=valstack(itop)

itop=itop-1

end subroutine popval

!**********************************************************************

SUBROUTINE GETPARAM_DC(sym,var)  ! Find double precision complex value var

                                 ! corresponding to symbol sym

character(len=*) :: sym

character(len=len_trim(sym)) :: usym

complex(kc8) :: var

ierr=0

sym=adjustl(sym)

if(is_letter(sym(1:1)).eqv..false. .or. len_trim(sym)>24) then

  ierr=11

  write(*,*) 'Error: symbol ',trim(sym),' has incorrect format'

  write(*,*)

  return

end if

ifind=0

usym=uppercase(sym)

do j=1,nparams

  if(trim(usym) == trim(params(j)%symbol)) then

    var=params(j)%value

    ifind=j

    exit

  end if

end do

if(ifind == 0) then          

  ierr=4

  write(*,*) 'Error: symbol ',trim(sym), ' not in symbol table'

  write(*,*) 

  return

end if

end subroutine getparam_dc

!**********************************************************************

SUBROUTINE GETPARAM_SC(sym,var)  ! Find single precision complex value var

                                 ! corresponding to symbol sym

character(len=*) :: sym

complex(kc4) :: var

complex(kc8) :: vard

call getparam_dc(sym,vard)

var=vard

end subroutine getparam_sc

!**********************************************************************

SUBROUTINE GETPARAM_DR(sym,var)  ! Find double precision real value var

                                 ! corresponding to symbol sym

character(len=*) :: sym

real(kr8) :: var

complex(kc8) :: vard

call getparam_dc(sym,vard)

var=real(vard,kr8)

end subroutine getparam_dr

!**********************************************************************

SUBROUTINE GETPARAM_SR(sym,var)  ! Find single precision real value var

                                 ! corresponding to symbol sym

character(len=*) :: sym

real(kr4) :: var

complex(kc8) :: vard

call getparam_dc(sym,vard)

var=real(vard)

end subroutine getparam_sr

!**********************************************************************

SUBROUTINE GETPARAM_DI(sym,ivar)  ! Find double precision integer value ivar

                                  ! corresponding to symbol sym

character(len=*) :: sym

integer(ki8) :: ivar

complex(kc8) :: vard

call getparam_dc(sym,vard)

ivar=nint(real(vard,kr8),ki8)

end subroutine getparam_di

!**********************************************************************

SUBROUTINE GETPARAM_SI(sym,ivar)  ! Find single precision integer value ivar

                                  ! corresponding to symbol sym

character(len=*) :: sym

integer(ki4) :: ivar

complex(kc8) :: vard

call getparam_dc(sym,vard)

ivar=nint(real(vard,kr8),ki4)

end subroutine getparam_si

!**********************************************************************

SUBROUTINE EVALEQN(eqn)  ! Evaluate an equation

character(len=*) :: eqn

character(len=len(eqn)) :: args(2)

complex(kc8) :: val

call parse(eqn,'=',args,nargs)   ! Seperate right- and left-hand-sides

call defparam(adjustl(args(1)),args(2)) ! Evaluate right-hand-side and

                                        ! assign to symbol on the

                                        ! left-hand-side.

end subroutine evaleqn

!**********************************************************************

SUBROUTINE LISTVAR      ! List all variables and their values

write(*,'(/a)') ' VARIABLE LIST:'

if(nparams == 0) then            ! Initialize symbol table

  params(1)%symbol='PI'

  params(1)%value=(3.14159265358979_kr8,0.0_kr8)

  params(2)%symbol='I'

  params(2)%value=(0.0_kr8,1.0_kr8)

  nparams=2

end if

do i=1,nparams

  write(*,*) trim(params(i)%symbol),' = ',params(i)%value

end do

end subroutine listvar

!**********************************************************************

end module evaluate