!*************************************************************************************!
! TITLE: CFVR4907 - PS 92: Seismic Analysis (design spectrum) Z-axis
! SUBTITLE: Seismic Analysis on a 1 DOF beam according to PS 92
!
! DESCRIPTION: This example compares spectrum accelerations given by CivilFEM to accelerations
! DESCRIPTION: obtained when solving on a 1 DOF beam and to the correct PS 92 spectum accelerations.
! DESCRIPTION:
! DESCRIPTION: The model is 1 element beam with only 1 DOF, which direction the CivilFEM spectrum
! DESCRIPTION: is set (the rest of the DOF are constrained). A different mass is applied at the
! DESCRIPTION: free end of the beam for the beam period to be equal to the 20 CivilFEM spectrum points.
! DESCRIPTION: After solving, the mass acceleration is obtained and compared to the CivilFEM spectrum
! DESCRIPTION: acceleration and to the correct value. Therefore, this example checks the following results:
! DESCRIPTION:
! DESCRIPTION: - CivilFEM spectrum is defined correctly by comparing to the correct value
! DESCRIPTION:
! DESCRIPTION: - CivilFEM spectrum is applied correctly and response acceleration is equal to spectrum acceleration
!
! ELEMENT TYPE: BEAM4, MASS21
! MODULES:
! UNITS: User
! KEYWORD1: Seismic
! KEYWORD2: PS 92
!
!*************************************************************************************!

   FINISH
  ~CFCLEAR,,1
  NomFile='CFVR4907'
  /TITLE, %NomFile%, PS 92: Seismic Analysis (design spectrum) Z-axis

! ----------------------------------------------------------------------
! Model definition
! ----------------------------------------------------------------------
! CivilFEM Setup: Code & Units
~UNITS,,LENG,CM
~UNITS,,TIME,S
~UNITS,,FORC,KP
~CODESEL,,,,,PS92
/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,,

! 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,AR20,EX,1         		          ! Ex   : Elastic Modulus
  ~CFGET,AR22,SECTION,1,MECHPROP,A,,2     ! Area : Bar area
  ~CFGET,AR23,SECTION,1,MECHPROP,IZZ,,2   ! Izz  : Moment of inertia
  ~CFGET,AR24,SECTION,1,MECHPROP,IYY,,2   ! Iyy  : Moment of inertia

! Spectrum definition
  ~DEFSPEC,ALL,0.1*981,S2,DESIGN,1,5

  pi=3.141592654

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

! Matrix 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,ARG1,1,20
   /PREP7
   ~CFGET,AR26,SEISM,,SPEC,TH,,ARG1,2   ! Period
  ! Labels
   LABEL(ARG1,1) = 'Sa(%ARG1%)z'
   LABEL(20+ARG1,1) = 'Sa(%ARG1%)z'

   AR25=12*AR20*AR24/L**3      ! Bending rigidity around Y axis
   AR27=AR25*(AR26/(2*pi))**2  ! Mass applied at end
   R,2,0,0,AR27
   TYPE,2
   MAT, 1
   REAL,2
   EN,2,10
   FINISH
   /SOLU

 ! Solution
 ! ---------------------------------------------------------------------
 ! Displacements
   D,1,,,,,,ALL
   D,10,,,,,UX,UY,,ROTX,ROTY,ROTZ
   ~MODLSOL,1
   ~CMBMOD,NONE,NONE,HORIZONT,0,0.0001
 ! Query results
   *GET,ARG5,NODE,10,U,Z
 ! Correct values
   VALUE(20-ARG1+1,1)=ARG5*AR25/AR27 ! Acceleration
 ! Obtained values
   ~CFGET,VALUE(20-ARG1+1,2),SEISM,,SPEC,SDH,,ARG1,2	 ! Y components of design horizontal spectrum
   ~CFGET,VALUE(40-ARG1+1,2),SEISM,,SPEC,SDH,,ARG1,2	 ! Y components of design horizontal spectrum
   /PREP7
   ! CivilFEM gives the spectrum Normalized by gravity.
   ! To compare it must be changed to acceleration user units.
   VALUE(20-ARG1+1,2) = VALUE(20-ARG1+1,2)*981
   VALUE(40-ARG1+1,2) = VALUE(40-ARG1+1,2)*981
   EDELE,2
*ENDDO

! Correct values
! ----------------------------------------------------------------------
 VALUE(40,1)= 13.01916
 VALUE(39,1)= 14.64925
 VALUE(38,1)= 16.63230
 VALUE(37,1)= 19.08261
 VALUE(36,1)= 22.16619
 VALUE(35,1)= 26.13204
 VALUE(34,1)= 31.36869
 VALUE(33,1)= 38.51187
 VALUE(32,1)= 48.66540
 VALUE(31,1)= 63.89720
 VALUE(30,1)= 68.23349
 VALUE(29,1)= 73.38913
 VALUE(28,1)= 79.64364
 VALUE(27,1)= 87.42696
 VALUE(26,1)= 97.44034
 VALUE(25,1)= 110.91976
 VALUE(24,1)= 130.29248
 VALUE(23,1)= 161.16594
 VALUE(22,1)= 220.72500
 VALUE(21,1)= 220.72500

! Warning and error tolerances
! ----------------------------------------------------------------------
   *DO,II,1,NComp
     TOLER(II, 1)= 0.1 $ TOLER(II, 2)= 0.1
   *ENDDO

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