dyn_matrix_test00prepnc.f90

Path: dynamics/dyn_matrix_test00prepnc.f90
Last Update: Thu Jun 21 12:23:59 JST 2007

dyn_matrix モジュールのテストプログラム用 NetCDF ファイル作成プログラム

Generate NetCDF file for test program for "dyn_matrix"

Authors:Yasuhiro MORIKAWA
Version:$Id: dyn_matrix_test00prepnc.f90,v 1.1.1.1 2007/06/21 03:23:59 morikawa Exp $
Tag Name:$Name: dcpam4-20070731 $
Copyright:Copyright (C) GFD Dennou Club, 2007. All rights reserved.

License::

Methods

Included Modules

dyn_spectral constants dc_args dc_types dc_string gt4_history

Public Instance methods

Main Program :

[Source]

program dyn_matrix_test00prepnc
  use dyn_spectral, only: DYNSP, Create, Close, PutLine, GetSpectralCoeff, GetDiffCoeff
  use constants, only: CONST, Create
  use dc_args, only: ARGS, Open, Debug, Help, Strict, Close
  use dc_types, only: DP, STRING
  use dc_string, only: StoA, toChar
  use gt4_history, only: GT_HISTORY, HistoryGet, HistoryCopy, HistoryCreate, HistoryAddVariable, HistoryPut, HistoryClose
  implicit none
  type(ARGS):: arg
  type(CONST):: c
  type(DYNSP):: dyn_sp00
  real(DP):: DelTime
  real(DP), allocatable:: wz_rn (:,:)
  integer, allocatable:: nmo(:,:,:)
  real(DP), allocatable:: wz_DiffVorDiv (:,:)
  real(DP), allocatable:: wz_DiffTherm (:,:)
  integer:: m, l, k, kk
  logical:: err
  character(STRING):: range_mes
  character(*), parameter:: output_file00 = 'dyn_matrix_test00.nc'
continue

  call Open(arg)
  call Debug(arg) 
   call Help(arg) 
   call Strict(arg)
  call Close(arg)

  call Create(c)

  DelTime = 180.0_DP

  call Create( dyn_sp = dyn_sp00, nmax = 10, imax = 32, jmax = 16, kmax = 2, PI = c % PI, RPlanet = c % RPlanet, Grav = c % Grav, Omega = c % Omega, Cp = c % Cp, RAir = c % RAir, VisOrder = 8, EFoldTime = 10800.0_DP, DelTime = DelTime ) ! (in)
  if ( allocated(nmo) ) deallocate(nmo)
  allocate( nmo(1:2, 0:10, 0:10) )
  if ( allocated(wz_rn) ) deallocate(wz_rn)
  allocate( wz_rn(11**2, 0:1) )
  call GetSpectralCoeff( dyn_sp = dyn_sp00, nmo = nmo, wz_rn = wz_rn ) ! (out)

  if ( allocated(wz_DiffVorDiv) ) deallocate(wz_DiffVorDiv)
  allocate( wz_DiffVorDiv(11**2, 0:1) )
  if ( allocated(wz_DiffTherm) ) deallocate(wz_DiffTherm)
  allocate( wz_DiffTherm(11**2, 0:1) )
  call GetDiffCoeff( dyn_sp = dyn_sp00, wz_DiffVorDiv = wz_DiffVorDiv, wz_DiffTherm = wz_DiffTherm ) ! (out)

  call PutLine( dyn_sp = dyn_sp00 )

  call HistoryCreate( file=output_file00, title='Spectral data for dyn_matrix_test', source='Spectral data generation program for dyn_matrix_test', institution='dcpam project / GFD Dennou Club', dims=StoA('m','l'), dimsizes=(/11, 11/), longnames=StoA('zonal wavenumber','meridional wavenumber'), units=StoA('1','1'), xtypes=StoA('int', 'int') )
  call HistoryPut('m', (/0.0,1.0,2.0,3.0,4.0,5.0,6.0,7.0,8.0,9.0,10.0/))
  call HistoryPut('l', (/0.0,1.0,2.0,3.0,4.0,5.0,6.0,7.0,8.0,9.0,10.0/))

  call HistoryAddVariable( varname = 'nmo1', dims = StoA('m', 'l'), longname = 'Spectral subscript expression 1', units='1', xtype='int')

  call HistoryAddVariable( varname = 'nmo2', dims = StoA('m', 'l'), longname = 'Spectral subscript expression 2', units='1', xtype='int' )

  call HistoryAddVariable( varname = 'DiffVorDiv1', dims = StoA('m', 'l'), longname = 'Coefficient of horizontal momentum diffusion 1 ', units='1/s', xtype='double' )

  call HistoryAddVariable( varname = 'DiffVorDiv2', dims = StoA('m', 'l'), longname = 'Coefficient of horizontal momentum diffusion 2 ', units='1/s', xtype='double' )

  call HistoryAddVariable( varname = 'DiffTherm1', dims = StoA('m', 'l'), longname = 'Coefficient of horizontal thermal and water diffusion 1', units='1/s', xtype='double' )

  call HistoryAddVariable( varname = 'DiffTherm2', dims = StoA('m', 'l'), longname = 'Coefficient of horizontal thermal and water diffusion 2', units='1/s', xtype='double' )

  call HistoryPut('nmo1', real(nmo(1,:,:)))
  call HistoryPut('nmo2', real(nmo(2,:,:)))

!!$  do l = 0, 10
!!$    do m = 0, 10
!!$      write(*,*) 'nmo(1,', m, ',', l, ')= ', nmo(1,m,l)
!!$    end do
!!$  end do
!!$
!!$  do l = 0, 10
!!$    do m = 0, 10
!!$      write(*,*) 'nmo(2,', m, ',', l, ')= ', nmo(2,m,l)
!!$    end do
!!$  end do

  do l = 0, 10
    do m = 0, 10
      call HistoryPut('DiffVorDiv1', dyn_sp00 % wz_DiffVorDiv(nmo(1,m,l),0), range='m='//trim(toChar(m))//',l='//trim(toChar(l)))
      if (m+l > 10) cycle
!!$      write(*,*) 'wz_DiffVorDiv(', 'nmo(1,', m, ',', l, '),0)= ', &
!!$        & dyn_sp00 % wz_DiffVorDiv(nmo(1,m,l),0)
    end do
  end do

  do l = 0, 10
    do m = 0, 10
      call HistoryPut('DiffVorDiv2', dyn_sp00 % wz_DiffVorDiv(nmo(2,m,l),0), range='m='//trim(toChar(m))//',l='//trim(toChar(l)))
      if (m+l > 10) cycle
!!$      write(*,*) 'wz_DiffVorDiv(', 'nmo(2,', m, ',', l, '),0)= ', &
!!$        & dyn_sp00 % wz_DiffVorDiv(nmo(2,m,l),0)
    end do
  end do

  do l = 0, 10
    do m = 0, 10
      call HistoryPut('DiffTherm1', dyn_sp00 % wz_DiffTherm(nmo(1,m,l),0), range='m='//trim(toChar(m))//',l='//trim(toChar(l)))
      if (m+l > 10) cycle
!!$      write(*,*) 'wz_DiffTherm(', 'nmo(1,', m, ',', l, '),0)= ', &
!!$        & dyn_sp00 % wz_DiffTherm(nmo(1,m,l),0)
    end do
  end do

  do l = 0, 10
    do m = 0, 10
      call HistoryPut('DiffTherm2', dyn_sp00 % wz_DiffTherm(nmo(2,m,l),0), range='m='//trim(toChar(m))//',l='//trim(toChar(l)))
      if (m+l > 10) cycle
!!$      write(*,*) 'wz_DiffTherm(', 'nmo(2,', m, ',', l, '),0)= ', &
!!$        & dyn_sp00 % wz_DiffTherm(nmo(2,m,l),0)
    end do
  end do

  call HistoryClose


end program dyn_matrix_test00prepnc

[Validate]