!*************************************************************************************!
! TITLE: CFVR4211 - GB50011-2010: Spectrum definition
! SUBTITLE: Spectrum definition according to GB 50011-2010
!
! DESCRIPTION: CivilFEM spectrum is defined and compared to the correct values.
!
! ELEMENT TYPE: NONE
! MODULES:
! UNITS: User
! KEYWORD1: GB 50011-2010
!
!*************************************************************************************!

   FINISH
  ~CFCLEAR,,1
  NomFile='CFVR4211'
  /TITLE, %NomFile%, GB50011-2010: Spectrum definition

! ----------------------------------------------------------------------
! Model definition
! ----------------------------------------------------------------------
! CivilFEM Setup: Code & Units
  ~UNITS,,LENG,CM
  ~UNITS,,TIME,S
  ~UNITS,,FORC,KP
  ~CODESEL,,,,,GB50011N
/PREP7
! Preprocessor
! ----------------------------------------------------------------------
! Materials: A-42
 ~CFMP,1,LIB,STEEL,EA,A42
! Modify density (Rho = 0)
 ~CFMP,1,USER
 ~CFMP,1,DatGen,RHO,,0
! Element Types
  ET,1,BEAM4   ! Type 1: 3D Beam
  ET,2,MASS21  ! Type 2: Mass
! Sections
 ~SSECLIB,1,1,1,1   !IPE 80
 ~BMSHPRO,1,BEAM,1,1,,,4,,0,,Steel

! Nodes
  L = 10              ! L    : Bar length
  N, 1
  N,10,L

! Elements
  TYPE,1
  MAT ,1
  REAL,1
  EN,1,1,10
  EPLOT

! Parameters
  *GET,_Ex,EX,1                         ! Ex   : Elastic Modulus
  ~CFGET,_A,SECTION,1,MECHPROP,A,,2    ! Area : Bar area
  ~CFGET,_Izz,SECTION,1,MECHPROP,IZZ,,2  ! Izz  : Moment of inertia
  ~CFGET,_Iyy,SECTION,1,MECHPROP,IYY,,2  ! Iyy  : Moment of inertia

! Seismic spectrum definition
 ~DEFSPEC,ALL,0.05,RARE,I0,3,4e-002,0.7,,

  pi=3.141592654

!  ----------------------------------------------------------------------
! DATA CHECK
! ----------------------------------------------------------------------
! Data comparision number
  NComp = 40
  NComp_ch = 0

! Marix dim.
  *DIM,LABEL,CHAR,Ncomp,1
  *DIM,LABEL_CH,CHAR,Ncomp_ch,1
  *DIM,VALUE,,Ncomp,3
  *DIM,VALUE_CH,CHAR,Ncomp_ch,3
  *DIM,TOLER,,Ncomp,2

*DO,II,1,20
   /PREP7
   ~CFGET,_T,SEISM,,SPEC,TH,,II,1   ! Period
  ! Labels
   LABEL(II,1) = 'Sa(%II%)x'
   LABEL(20+II,1) = 'Sa(%II%)x'

   _K=_Ex*_A/L
   _M=_K*(_T/(2*pi))**2  ! Mass applied at end
   R,2,_M
   TYPE,2
   MAT, 1
   REAL,2
   EN,2,10
   FINISH
   /SOLU

 ! Displacements
   D,1,,,,,,ALL
   D,10,,,,,,UY,UZ,ROTX,ROTY,ROTZ
   ~MODLSOL,1
   ~CMBMOD,HORIZONT,NONE,NONE,0,0.0001
 ! Query results
   *GET,_Ydisp,NODE,10,U,X
 ! Correct values
   VALUE(20-II+1,1)=_Ydisp*_K/_M ! Aceleration
 ! Obtained values
   ~CFGET,VALUE(20-II+1,2),SEISM,,SPEC,SDH,,II,1  ! Y components of frequent horizontal spectrum
   ~CFGET,VALUE(40-II+1,2),SEISM,,SPEC,SDH,,II,1  ! Y components of frequent horizontal spectrum
   
   ! CivilFEM gives the spectrum Normalized by gravity.
   ! To compare it must be changed to acceleration user units.
   VALUE(20-II+1,2) = VALUE(20-II+1,2)*981
   VALUE(40-II+1,2) = VALUE(40-II+1,2)*981
   /PREP7
   EDELE,2
*ENDDO

! Correct values
! ----------------------------------------------------------------------
  VALUE(21,1) = 140.6209
  VALUE(22,1) = 157.6358
  VALUE(23,1) = 174.6507
  VALUE(24,1) = 191.6656
  VALUE(25,1) = 208.6805
  VALUE(26,1) = 225.6954
  VALUE(27,1) = 242.7103
  VALUE(28,1) = 259.7252
  VALUE(29,1) = 276.7401
  VALUE(30,1) = 293.755
  VALUE(31,1) = 293.755
  VALUE(32,1) = 183.7439645
  VALUE(33,1) = 134.8937338
  VALUE(34,1) = 107.0879388
  VALUE(35,1) = 89.06180838
  VALUE(36,1) = 76.39279111
  VALUE(37,1) = 66.98358301
  VALUE(38,1) = 57.9628421
  VALUE(39,1) = 48.94210119
  VALUE(40,1) = 39.92136028
                        
! Warning and error tolerances
! ----------------------------------------------------------------------
   *DO,II,1,40
     TOLER(II, 1)= 1.0 $ TOLER(II, 2)= 1.0
   *ENDDO

! ----------------------------------------------------------------------
! Results Comparison
!-----------------------------------------------------------------------
  COMPARA.MAC