! copyright 2005, J. E. Akin, all rights reserved. MODULE PRECISION_MODULE IMPLICIT NONE ! Use hardware defaults for single and double precision INTEGER, PARAMETER :: SP = KIND (1.0) ! Single precision INTEGER, PARAMETER :: DP = KIND (1.d0) ! Double precision ! Note: Use SELECTED_INT_KIND or SELECTED_REAL_KIND for user ! defined portable precision. However, user choice may not be ! on all hardware so check it with functionality below. CONTAINS SUBROUTINE CHECK_HARDWARE_PRECISION ! Verify that the hardware supports the request IF ( KIND (SP) < 0 ) THEN PRINT *, 'ERROR: SP PRECISION NOT AVAILABLE' STOP 'SET SP = KIND (1.0) IN MODULE PRECISION_MODULE' END IF IF ( KIND (DP) < 0 ) THEN PRINT *, 'ERROR: DP PRECISION NOT AVAILABLE' STOP 'SET DP = KIND (1.d0) IN MODULE PRECISION_MODULE' END IF END SUBROUTINE CHECK_HARDWARE_PRECISION END MODULE PRECISION_MODULE Module System_Constants Use Precision_Module IMPLICIT NONE Integer, parameter :: FOURIER_NEG = -1 Integer :: & DATA_SET, DEG_HOMO, EIGEN_POST, EIGEN_SCP, EIGEN_SHOW, & EXACT_CASE, EXACT_FL, EXACT_FX, & EXACT_INT, EXACT_REALS, EXAMPLE, GET_DAT, & HISTORY_DOF, HISTORY_NODE, HISTORY_PARM, INC_SAVE, & IN_RHS, IP_TEST, IS_MIX, IS_SEG, IS_STD, & IS_USR, I_BUG, I_SAY, LEM_WRT, LINE_IN, LINE_WARN, & LIST_LR, & LP_TEST, L_B_N, L_CLASS, L_CONT, L_HOMO, L_REACT, & L_SHAPE, L_S_TOT, MAT_FLO, MAX_ACT, MAX_DAT, MAX_FACES,& MAX_NP, & MAX_TYP, MISC_FL, MISC_FX, MODE, M_B_N, M_COEFF, & M_SHOW, NEEDS, NEIGH_L, NEIGH_N, NEIGH_P, NEXT_SAVE,& NOD_PER_EL, NOT_SYM, NO_DIST, & NPT_WRT, NUL_COL, N_ANAL, & N_BS_FIX, N_BS_FLO, N_BUG, N_CEQ, & N_COEFF, N_CORNER, N_CRD, N_D_FRE, N_ELEMS, N_EL_FRE, & N_D_FLUX, N_GEOM, N_G_DOF, N_G_FLUX, N_HOMO, N_ITER, & N_LP_FIX, N_LP_FLO, N_L_DEL, N_L_TYPE, N_MAT, N_MIXED, & N_MX_FIX, N_MX_FLO, N_NP_FIX, N_NP_FLO, N_PARM, N_PATCH, & N_PRT, N_QP, N_QP_C, N_QP_R, N_REACT, N_R_B, & N_F_SEG, N_SPACE, N_STORE, & N_FILE1, N_FILE2, N_FILE3, N_FILE4, N_FILE5, N_T_SUM, & N_WARN, N_2_DER, N_3_DER, PRT_STEP, RESTART_STEP, & SCP_DEG, SCP_FIT, SCP_INC, SCP_PAD, SCP_QP, SCP_RECL, & SWEEPS, THIS_EL, THIS_ITER, & THIS_LT, THIS_STEP, TIME_METHOD, & T_SETS, T_STEPS, U_DOF, U_EL_ADJ, U_FLUX, U_HIST, & U_INTGR, U_LOG, U_NODE, U_ON_B, & U_PLT1, U_PLT2, U_PLT3, U_PLT4, & U_PLT5, U_PLT6, U_PLT7, U_PLT8, U_PLT9, U_PT_ADJ, & U_SCPR, U_TIME, U_2_DER, U_3_DER, U2_PLT1, U2_PLT2, & U2_PLT3, U2_PLT4, U_GRAD !b convert all U_GRAD to U_FLUX '123456789012', 123456789012 Integer :: & ! property counts el_int, el_real, integers, mat_real, mixed_int, mixed_real, & pt_int, pt_real, reals, seg_int, seg_real ! variable dof per node: DOF_PT (0:MAX_NP), DOF_PT_SUM (-1:MAX_NP) INTEGER, ALLOCATABLE :: DOF_PT (:) ! number of dof at each node INTEGER, ALLOCATABLE :: DOF_PT_SUM (:) ! total eqs, including this node ! EQ_IJ = EQ NUMBER: DOF_PT_SUM (I-1) + J, J <= DOF_PT(I), I,J > 0 CHARACTER(len=12) :: RECTANGULAR (3) = (/' X-Coord, ', & ' Y-Coord, ', & ' Z-Coord, ' /) CHARACTER(len=12) :: CYLINDRICAL (3) = (/' Radius r, ', & ' Axial z, ', & ' Angle xz_r,' /) CHARACTER(len=12) :: POLAR_SPACE (4) = (/' Radius r, ', & ' Angle_to_X,', & ' Angle_to_Y,', & ' Angle_to_Z,' /) CHARACTER(len=12) :: SPHERICAL_S (3) = (/' Radius r, ', & ' Ang_to_xz, ', & ' Ang_to_z, ' /) CHARACTER(len=10) :: RECT_SHORT (3) = (/' X-Coord, ', & ' Y-Coord, ', & ' Z-Coord, ' /) CHARACTER(len=10) :: CYL_SHORT (3) = (/' Radius r,', & ' Axial z, ', & ' Ang xz_r,' /) CHARACTER(len=10) :: SPHER_SHORT (3) = (/' Radius r,', & ' Ang_xz, ', & ' Ang_to_z,' /) CHARACTER(len=10) :: POLAR_SHORT (4) = (/' Radius r,', & ' Ang_to_X,', & ' Ang_to_Y,', & ' Ang_to_Z,' /) CHARACTER(len=12), ALLOCATABLE :: XYZ_NAME (:), DOF_NAME (:), & EXACT_NAME (:), FLUX_NAME (:), NO_NAME (:), PROP_NAME (:), & DERIV_2_NAME (:), DERIV_3_NAME (:), AVE_NAME (:), INPUT_NAME (:) CHARACTER(len=10), ALLOCATABLE :: XYZ_SHORT (:), DOF_SHORT (:), & EXACT_SHORT (:), FLUX_SHORT (:), NO_SHORT (:), PROP_SHORT (:),& DERIV_2_SHORT (:), DERIV_3_SHORT (:), AVE_SHORT (:) CHARACTER(len=70) :: TITLE CHARACTER(len=3) :: STEP_STR ! step counted as 3 char string REAL(DP), PARAMETER :: TWO_PI = 6.28318530717958623d0 REAL(DP), PARAMETER :: ONE_PI = TWO_PI * 0.5d0 REAL(DP), PARAMETER :: PI_SQ = ONE_PI * ONE_PI REAL(DP), PARAMETER :: PI = ONE_PI REAL(DP) :: AREA_THICK, CUT_OFF, MAX_DERIVATIVE REAL(DP) :: PERCENT_ERR_MAX, RELAXATION, SCALAR_SOURCE ! Iterative load control REAL(DP) :: LOAD_STEPS, LOAD_INC ! Eigenvalue and/or timestepping !b REAL(DP) :: DENSITY !b ! Time stepping: transient & dynamic REAL(DP) :: TIME, TIME_STEP, START_TIME, TS_B, TS_A REAL(DP) :: INITIAL_VALUE, TRAPEZOIDAL_COEFF REAL(DP) :: FREQ_COS, FREQ_SIN REAL(DP) :: STEP_UP_TIME, STEP_DOWN_TIME REAL(DP) :: RAMP_UP_TIME, RAMP_DOWN_TIME REAL(DP) :: RAMP_UP_TIMES (2), RAMP_DOWN_TIMES (2) REAL(DP) :: UP_DOWN_TIMES (2) ! Dynamic REAL(DP) :: HUGHES_ALPHA, NEWMARK_BETA, NEWMARK_GAMMA REAL(DP) :: WILSON_THETA, MASS_DAMPING, STIF_DAMPING ! Normal flux on a boundary segment: K_n * U,n + Q_NORMAL_SEG = 0 REAL(DP) :: Q_NORMAL_SEG, FLUX_THICK, THICKNESS_SEG ! Mixed or Robin boundary condition: ! K_n * U,n + ROBIN_1_SEG * U + ROBIN_2_SEG = 0, or special ! K_n * U,n + CONVECT_COEF * ( U - CONVECT_TEMP) = 0 ! over optional line thickness CONVECT_THICK, or ROBIN_THICK REAL(DP) :: CONVECT_COEF, CONVECT_TEMP, CONVECT_THICK REAL(DP) :: ROBIN_1_SEG, ROBIN_2_SEG, ROBIN_THICK REAL(DP) :: AVE_H_WT ! geometric bounds REAL(DP), ALLOCATABLE :: MAX_XYZ (:) ! max mesh bounds REAL(DP), ALLOCATABLE :: MIN_XYZ (:) ! min mesh bounds ! dynamic analysis REAL(DP), ALLOCATABLE :: DD_DOT (:) ! system velocities REAL(DP), ALLOCATABLE :: DD_2DOT (:) ! system accelerations ! eigenvector post processing REAL(DP) :: THIS_SCALE LOGICAL :: ANISOTROPIC = .FALSE. LOGICAL :: AVERAGE_MASS = .FALSE. LOGICAL :: AXISYMMETRIC = .FALSE. LOGICAL :: BAR_CHART = .FALSE. LOGICAL :: BUCKLING = .FALSE. LOGICAL :: CONVECTION = .FALSE. LOGICAL :: CONVECT_VARY = .FALSE. LOGICAL :: CONSTANT_J = .FALSE. LOGICAL :: DATA_CHECK_ONLY = .FALSE. LOGICAL :: DEBUG_ADAPT = .FALSE. LOGICAL :: DEBUG_ASSEMBLY = .FALSE. LOGICAL :: DEBUG_B = .FALSE. LOGICAL :: DEBUG_E = .FALSE. LOGICAL :: DEBUG_EL_COL = .FALSE. LOGICAL :: DEBUG_EL_SQ = .FALSE. LOGICAL :: DEBUG_EL_POST = .FALSE. LOGICAL :: DEBUG_EL_TYPE = .FALSE. LOGICAL :: DEBUG_ERROR = .FALSE. LOGICAL :: DEBUG_EXACT = .FALSE. LOGICAL :: DEBUG_INERTIA = .FALSE. LOGICAL :: DEBUG_INCLUDE = .FALSE. LOGICAL :: DEBUG_MIX_SQ = .FALSE. LOGICAL :: DEBUG_POST_EL = .FALSE. LOGICAL :: DEBUG_PROPERTY = .FALSE. LOGICAL :: DEBUG_SCP = .FALSE. LOGICAL :: DEBUG_SEG_COL = .FALSE. LOGICAL :: DEBUG_UNITS = .FALSE. LOGICAL :: DIAGONAL_MASS = .FALSE. LOGICAL :: DOF_VARY = .FALSE. LOGICAL :: DYNAMIC = .FALSE. LOGICAL :: ECHO_BC = .TRUE. LOGICAL :: ECHO_EL = .TRUE. LOGICAL :: ECHO_FLUX = .TRUE. LOGICAL :: ECHO_INPUT = .TRUE. LOGICAL :: ECHO_MIXED = .TRUE. LOGICAL :: ECHO_PROP = .TRUE. LOGICAL :: ECHO_PTS = .TRUE. LOGICAL :: ECHO_RHS = .TRUE. LOGICAL :: ECHO_START = .TRUE. LOGICAL :: EIGEN = .FALSE. LOGICAL :: EIGEN_ROOT = .FALSE. LOGICAL :: EL_DAMPING = .FALSE. LOGICAL :: FLUX_NORMAL = .FALSE. LOGICAL :: FORCE_BODY = .FALSE. LOGICAL :: FOURIER_LAW = .FALSE. LOGICAL :: FREQ_CONTROL = .FALSE. LOGICAL :: FROM_REST = .FALSE. LOGICAL :: F90 = .FALSE. LOGICAL :: F95 = .FALSE. LOGICAL :: GEOM_ITER = .FALSE. LOGICAL :: GEOMETRIC_K = .FALSE. LOGICAL :: GET_DELETED_EL = .FALSE. LOGICAL :: GET_2ND_DERIV = .FALSE. LOGICAL :: GLOBAL_PROPERTY = .FALSE. LOGICAL :: GRAD_BASE_ERROR = .FALSE. LOGICAL :: HAS_BODY_FORCE = .FALSE. LOGICAL :: HEAVYSIDE = .FALSE. LOGICAL :: HUGHES_METHOD = .FALSE. LOGICAL :: INCREMENTAL = .FALSE. LOGICAL :: INITIAL_FUNCTION = .FALSE. LOGICAL :: INPUT_ONLY = .FALSE. LOGICAL :: IS_ELEMENT = .FALSE. LOGICAL :: IS_FLUX_SEG = .FALSE. LOGICAL :: IS_MIXED_BC = .FALSE. LOGICAL :: IS_NULL = .FALSE. LOGICAL :: JACOBI = .FALSE. LOGICAL :: LIST_EL_TO_EL = .FALSE. LOGICAL :: LIST_QP_FLUX = .FALSE. LOGICAL :: NEWMARK_METHOD = .TRUE. LOGICAL :: NODAL_TAU = .FALSE. LOGICAL :: NO_ERROR_EST = .FALSE. LOGICAL :: NO_PRINTING = .FALSE. LOGICAL :: NO_PROPERTIES = .FALSE. LOGICAL :: NO_SCP_AVE = .FALSE. LOGICAL :: NULL_SCP_AVE = .FALSE. LOGICAL :: PLANE_STRAIN = .FALSE. LOGICAL :: POLAR_COORD = .FALSE. LOGICAL :: POST_EL = .FALSE. LOGICAL :: POST_MIXED = .FALSE. LOGICAL :: PROPERTY_FIRST = .FALSE. LOGICAL :: PT_EL_LIST = .FALSE. LOGICAL :: PT_EL_SUMS = .FALSE. LOGICAL :: PT_MASSES = .FALSE. LOGICAL :: PT_DAMPING = .FALSE. LOGICAL :: RAMP_CONTROL = .FALSE. LOGICAL :: RAMP_ITER = .FALSE. LOGICAL :: SAVE_NEW_MESH = .FALSE. LOGICAL :: SAVE_1248 = .FALSE. LOGICAL :: SCALE_FREQ = .FALSE. LOGICAL :: SCALE_LATE = .FALSE. LOGICAL :: SCALE_LOADS = .FALSE. LOGICAL :: SCALE_RAMP = .FALSE. LOGICAL :: SCALE_STEP = .FALSE. LOGICAL :: SCALE_ZERO = .FALSE. LOGICAL :: SCP_CENTER_ONLY = .FALSE. LOGICAL :: SCP_DOF = .FALSE. LOGICAL :: SCP_E_AT_NP = .FALSE. LOGICAL :: SCP_E_AT_QP = .FALSE. LOGICAL :: SCP_E_SMOOTH = .FALSE. LOGICAL :: SCP_NEIGH_EL = .TRUE. LOGICAL :: SCP_NEIGH_FACE = .FALSE. LOGICAL :: SCP_NEIGH_PT = .FALSE. LOGICAL :: SCP_ONLY_ONCE = .TRUE. LOGICAL :: SCP_2ND_DERIV = .FALSE. LOGICAL :: SCP_3RD_DERIV = .FALSE. LOGICAL :: SCP_REC_NUMB_ALLOC = .FALSE. LOGICAL :: SCP_RECORDS_SAVED = .FALSE. LOGICAL :: SKIP_ERROR = .FALSE. LOGICAL :: SPACE_TIME = .FALSE. LOGICAL :: STATIC = .TRUE. LOGICAL :: SUPG = .FALSE. LOGICAL :: SYMMETRIC = .TRUE. LOGICAL :: S_T_CONTINUOUS = .FALSE. LOGICAL :: TAU_BOX = .FALSE. LOGICAL :: TAU_GEOM = .FALSE. LOGICAL :: TAU_PEC = .FALSE. LOGICAL :: TAU_QP = .FALSE. LOGICAL :: TAU_SHOW = .FALSE. LOGICAL :: TAU_S1 = .FALSE. LOGICAL :: TAU_RGN = .FALSE. LOGICAL :: TAU_UGN = .FALSE. LOGICAL :: TAU_VOL = .FALSE. LOGICAL :: TIME_SLAB = .FALSE. LOGICAL :: TOUCH = .FALSE. LOGICAL :: TRANSIENT = .FALSE. LOGICAL :: TYPES_ALLOCATED = .FALSE. LOGICAL :: UPWIND = .FALSE. LOGICAL :: USE_EXACT = .FALSE. LOGICAL :: USE_EXACT_BC = .FALSE. LOGICAL :: USE_EXACT_FLUX = .FALSE. LOGICAL :: USE_EXACT_ROBIN = .FALSE. LOGICAL :: USE_EXACT_SOURCE = .FALSE. LOGICAL :: USER_LOGIC = .FALSE. LOGICAL :: VARY_E = .FALSE. LOGICAL :: WILSON_METHOD = .FALSE. ! System Arrays REAL(DP) :: BODY_FORCE (3) ! needed in keywords !b REAL(DP), ALLOCATABLE :: BODY_F (:) ! will be body !b REAL(DP), ALLOCATABLE :: X (:,:) !b integer, ALLOCATABLE :: I_BC (:), NODES (:,:), N_RES_EQ (:) ! Friendly elem type data controls INTEGER, ALLOCATABLE :: TYPE_NODES (:), TYPE_GAUSS (:), & TYPE_GEOM (:), TYPE_PARM (:), TYPE_SHAPE (:), & TYPE_FACES (:), TYPE_SONS (:), TYPE_USER (:) INTEGER :: RECORD_NUMBER INTEGER, ALLOCATABLE :: SCP_RECORD_NUMBER (:, :) CONTAINS ! ENCAPSULATED FUNCTIONALLY FUNCTION TO3STR (K) ! integer to 3 character string CHARACTER(LEN=3) :: TO3STR INTEGER, INTENT(IN) :: K ! for file names: K=2, TO3STR = '002' WRITE (TO3STR, '(I3.3)') MOD (K, 999) END FUNCTION TO3STR FUNCTION TO4STR (K) ! integer to 4 character string CHARACTER(LEN=4) :: TO4STR INTEGER, INTENT(IN) :: K ! for file names: K=2, TO4STR = '0002' WRITE (TO4STR, '(I4.4)') MOD (K, 9999) END FUNCTION TO4STR SUBROUTINE SET_TO_MIXED_OR_SEG_OR_STD (IE) ! .............................................................. ! SET ELEMENT GROUP, AS STORED IN NODES ARRAY ! .............................................................. INTEGER, INTENT(IN) :: IE ! Element number IS_MIX = 0 ; IS_SEG = 0 ; IS_STD = 0 ; IS_USR = 0 if ( IE <= N_ELEMS ) then IS_ELEMENT = .TRUE. ; IS_MIXED_BC = .FALSE. ; IS_FLUX_SEG = .FALSE. IS_STD = IE ; THIS_EL = IE else if ( IE > N_ELEMS .and. IE <= (N_ELEMS + N_MIXED) ) then IS_ELEMENT = .FALSE. ; IS_MIXED_BC = .TRUE. ; IS_FLUX_SEG = .FALSE. IS_MIX = IE - N_ELEMS ; THIS_EL = 0 else if ( IE > (N_ELEMS + N_MIXED) .and. & IE <= (N_ELEMS + N_MIXED + N_F_SEG)) then IS_ELEMENT = .FALSE. ; IS_MIXED_BC = .FALSE. ; IS_FLUX_SEG = .TRUE. IS_SEG = IE - (N_ELEMS + N_MIXED) ; THIS_EL = 0 else IS_USR = IE - (N_ELEMS + N_MIXED + N_F_SEG) print *,'WARNING: unknown user group for element ', IE N_WARN = N_WARN + 1 IS_ELEMENT = .FALSE. ; IS_MIXED_BC = .FALSE. ; IS_FLUX_SEG = .FALSE. end if END SUBROUTINE SET_TO_MIXED_OR_SEG_OR_STD SUBROUTINE SET_CONTROL_DEFAULTS ! * * * * * * * * * * * * * * * * * * * * * * * * * * ! Set control defaults to override with keyword inputs ! * * * * * * * * * * * * * * * * * * * * * * * * * * ! Characters TITLE = 'Untitled' ! Reals AREA_THICK = 1.d0 ! Thickness of 2-D domain BODY_FORCE = 0.d0 ! No body forces per unit volume CONVECT_COEF = 0.d0 ! Mixed convection data CONVECT_TEMP = 0.d0 ! Mixed convection data CONVECT_THICK = 1.d0 ! Mixed convection optional line thickness CUT_OFF = 1.d-5 ! Iteration control DENSITY = 0.d0 ! Default mass density !b FLUX_THICK = 1.d0 ! Flux segment thickness FREQ_COS = 0.d0 ! Scale loads at Cos (FREQ_COS*TIME) FREQ_SIN = 0.d0 ! Scale loads at Sin (FREQ_SIN*TIME) HUGHES_ALPHA = -0.05d0 ! Default time integration constant MASS_DAMPING = 0.d0 ! Multiplies mass matrix MAX_DERIVATIVE = 0.d0 ! Plotting scales NEWMARK_BETA = 0.25d0 ! Default time integration constant NEWMARK_GAMMA = 0.5d0 ! Default time integration constant PERCENT_ERR_MAX = 1.d0 ! SCP default energy norm accuracy Q_NORMAL_SEG = 0.d0 ! Segment normal flux RAMP_UP_TIME = 0.d0 ! Ramp from 0 to 1 interval RAMP_DOWN_TIME = 0.d0 ! Ramp from 1 to 0 interval RAMP_UP_TIMES = (/ 0.d0, 1.d0 /) ! Ramp from 0 to 1 interval RAMP_DOWN_TIMES = (/ 0.d0, 1.d0 /) ! Ramp from 1 to 0 interval RELAXATION = 1.25d0 ! Iteration control ROBIN_1_SEG = 0.d0 ! Mixed Robin data ROBIN_2_SEG = 0.d0 ! Mixed Robin data ROBIN_THICK = 1.d0 ! Mixed Robin optional line thickness AVE_H_WT = 1.d0 ! New mesh scale, R_ave or 0 = R_min SCALAR_SOURCE = 0.d0 ! Scalar source per unit volume START_TIME = 0.d0 ! Time integration STIF_DAMPING = 0.d0 ! Multiplies stiffness matrix THICKNESS_SEG = 1.d0 ! Flux segment thickness THIS_SCALE = 1.d0 ! Current eigenvector normalization scale TIME = 0.d0 ! Time integration TIME_STEP = 1.d0 ! Time integration TRAPEZOIDAL_COEFF = 0.5d0 ! Gen Trap integration, Crank-Nicolson TS_B = 0.5d0 ! Time integration, Crank-Nicolson TS_A = 1.d0 - TS_B ! Time integration UP_DOWN_TIMES = (/ 0.5d0, 1.d0 /) ! Ramp from 0 to 1 to 0 interval ! Integers DATA_SET = 0 ! Data file number for "example" application DEG_HOMO = 0 ! Default to non_homogeneous polynomial degrees EIGEN_POST = 0 ! Number of eigenvectors for element postprocessing EIGEN_SCP = 1 ! The mode to use for scp and/or averaging fluxes EIGEN_SHOW = 10 ! Number of eigen values & vectors to list EXACT_CASE = 0 ! Previously created exact source case EXACT_FL = 0 ! Number of exact solution real properties EXACT_FX = 0 ! Number of exact solution integer properties EXAMPLE = 0 ! Previously created source "example" GET_DAT = 5 ! Number of items to read in element type definition HISTORY_DOF = 0 ! DOF number for time history graph HISTORY_NODE = 0 ! Node number for time history graph of a parameter HISTORY_PARM = 0 ! Parameter number for time history graph at a node INC_SAVE = 0 ! Increment of save step number, if > 0 IN_RHS = 0 ! Initial column matrix flag. Read if = 1. IP_TEST = 0 ! Property existance flag, = 0 if none IS_MIX = 0 ! Mixed bc element number IS_SEG = 0 ! Flux segment element number IS_STD = 1 ! Standard element number IS_USR = 0 ! Element number in user defined element group I_BUG = 0 ! Turn on debug prints if > 0 I_SAY = 0 ! Number of user remark lines LEM_WRT = 0 ! Element's nodal parameter print flag, =0 omit print LINE_IN = 0 ! Last input line read LINE_WARN = 5 ! Maximum warnings for one input line LIST_LR = 0 ! Element reactions flag, 0=Omit, 1=Compute & List LP_TEST = 0 ! Element properties flag, 0=None, else will read L_B_N = 0 ! Maximum nodes on an element boundary segment L_CLASS = 0 ! Undefined element class 1=el, 2=mixed_bc 3=flux seg L_HOMO = 0 ! Element properties flag, = 0 if non-homogeneous L_REACT = 0 ! Element reactions flag, 0=Omit, 1=List L_REACT = 21 ! Element reaction data storage unit: > 0 use & list L_SHAPE = 1 ! Element shape, 1=line L_S_TOT = 0 ! Total number of elements and segments, >=N_ELEMS MAT_FLO = 0 ! Number of real material properties MAX_ACT = 0 ! Number of active constraint types, <=MAX_TYP MAX_DAT = 8 ! Number of items in an element type definition MAX_FACES = 8 ! Maximum number of element faces MAX_NP = 2 ! Number of nodal points in the system MAX_TYP = 3 ! Max number of dof in linear constraint MISC_FL = 0 ! Number of miscellaneous real properties MISC_FX = 0 ! Number of miscellaneous integer properties MODE = 0 ! Unsymmetric storage flag 0=False, 1=True M_B_N = 0 ! Maximum nodes on a mixed_bc boundary segment = L_B_N M_COEFF = 0 ! Number of terms in system mass matrix M_SHOW = 20 ! Max assembled sq matrix debug print size NEEDS = 0 ! Number of shared nodes to be an element neighbor NEIGH_L = 0 ! Maximum number of element neighbors at any elem NEIGH_N = 0 ! Maximum number of element neighbors at any node NOD_PER_EL = 2 ! Maximum number of nodes per element NEXT_SAVE = 1 ! Next save step in 1, 2, 4, 8, 16, ... NOT_SYM = 0 ! Number of terms in lower skyline NO_DIST = 1 ! No distances in node plots NPT_WRT = 1 ! Nodal parameter print flag. List by nodes if 1. NUL_COL = 0 ! Element column matrix flag, =0 if required N_ANAL = 0 ! Number of element analysis nodes <= NOD_PER_EL N_BS_FIX = 0 ! Number of boundary segment integer properties N_BS_FLO = 0 ! Number of boundary segment real properties N_BUG = 6 ! Debug printer unit N_CEQ = 0 ! Number of active constraint equations N_COEFF = 0 ! Number of terms in upper skyline N_CORNER = 0 ! Number of corner nodes, <= NOD_PER_EL N_CRD = 5 ! Input unit (card reader) N_D_FRE = 2 ! Number of degrees-of-freedom in the system N_ELEMS = 1 ! Number of elements in the system N_EL_FRE = 2 ! Maximum number of degrees-of-freedom in an elem N_D_FLUX = 0 ! Maximum number of flux entries, = L_B_N*N_G_DOF N_GEOM = 0 ! Number of geometry nodes, <= NOD_PER_EL N_G_DOF = 1 ! Number of parameters per node N_G_FLUX = 0 ! Number of generalized flux components, if N_F_SEG>0 N_HOMO = 0 ! Nodal properties flag, = 0 if non-homogeneous N_ITER = 1 ! Number of iterations (negative turns on debug) N_LP_FIX = 0 ! Number of integer proterties per element N_LP_FLO = 0 ! Number of real properties per element N_L_DEL = 0 ! Number of deleted elements N_L_TYPE = 1 ! Number of element types N_MAT = 0 ! Number of materials N_MIXED = 0 ! Number of segments with mixed boundary conditions N_MX_FIX = 0 ! Number of integer properties per mixed segment N_MX_FLO = 0 ! Number of real properties per mixed segment N_NP_FIX = 0 ! Number of integer properties per node N_NP_FLO = 0 ! Number of real properties per node N_PARM = 0 ! Parametric geometry dimensions (see L_SHAPE) N_PATCH = 0 ! Number of SCP patches, MAX_NP or N_ELEMS N_PRT = 6 ! Line printer unit N_QP = 0 ! Maximum number of quadrature points in an element N_QP_C = 0 ! Maximum number of quadrature points for element C N_QP_R = 0 ! Maximum number of reduced integration points N_REACT = 20 ! Nodal reaction data storage unit N_R_B = 1 ! Number of rows in element B matrix N_F_SEG = 0 ! Number of element boundary segments with flux N_SPACE = 1 ! Dimension of solution space N_STORE = 0 ! Initialize storage size N_FILE1 = 0 ! Flag for post solution calculation data, else 32 N_FILE2 = 0 ! Flag for post solution calculation data, else 33 N_FILE3 = 0 ! Scratch unit, else 34 N_FILE4 = 0 ! Scratch unit, else 35 N_FILE5 = 0 ! Scratch unit, else 36 N_T_SUM = 0 ! Number of scratch units active N_WARN = 0 ! Number of warning messages issued N_2_DER = 1 ! Number of (exact) second derivatives, 1, 3, or 6 N_3_DER = 1 ! Number of (exact) third derivatives, 1, 4, or 9 RESTART_STEP = 0 ! Restart solution from this step SCP_DEG = 0 ! Default to non_homogeneous element degree SCP_FIT = 0 ! Number of components fit in SCP, =N_R_B SCP_INC = 0 ! Raise SCP degree this amount above element degree SCP_PAD = 0 ! default record pad length SCP_RECL = 0 ! default record length SWEEPS = 15 ! Max sweeps in general Jacobi eigensolution THIS_EL = 0 ! current element number THIS_ITER = 0 ! current iteration number THIS_STEP = 0 ! current time step TIME_METHOD = 2 ! Crank-Nicolson transient , or Newmark dynamic T_SETS = 1 ! Number of time integration groups T_STEPS = 0 ! Number of time integration steps in this group U_DOF = 58 ! Unit to store dof time-history U_EL_ADJ = 55 ! Optional element adjacent elements save unit U_PT_ADJ = 56 ! Optional nodal adjacent elements save unit U_FLUX = 30 ! Elem flux from EL_SQ to POST_APP & SCP U_GRAD = 32 ! Elem gradient from EL_SQ to POST_APP & SCP !b U_GRAD = U_FLUX !b remove U_GRAD from old el includes U_HIST = 0 !b 41 ! Unit to append miscellaneous job logs U_INTGR = 31 ! Elem H integral from EL_SQ U_LOG = 0 !b 42 ! Unit to store a single job log U_NODE = 59 ! Unit to store nodal time history U_ON_B = 57 ! Optional boundary faces for plotters U_PLT1 = 22 ! Optional coord output for plotters U_PLT2 = 23 ! Optional topology output for plotters U_PLT3 = 24 ! Optional solution output for plotters U_PLT4 = 25 ! Optional elem gradient output for plotters U_PLT5 = 26 ! Optional node gradient output for plotters U_PLT6 = 27 ! Optional element error output for plotters U_PLT7 = 28 ! Optional pt ave error output for plotters U_PLT8 = 34 ! Optional exact solution output for plotters U_PLT9 = 35 ! Optional exact flux output for plotters U_SCPR = 29 ! Active SCP direct file unit U_TIME = 60 ! Unit for a time history graph U2_PLT1 = 51 ! Optional fe grad of fluxes output for plotters U2_PLT2 = 52 ! Optional exact grad of fluxes output for plotters U2_PLT3 = 53 ! Optional 2nd derivative data output for plotters U2_PLT4 = 54 ! Optional 2nd derivative data output for plotters U_2_DER = 37 ! Optional SCP 2nd derivative unit U_3_DER = 38 ! Optional SCP 3rd derivative unit ! Integer property input counts el_int = 0 ! Number of integer properties per element el_real = 0 ! Number of real properties per element integers = 0 ! Number of miscellaneous integer properties mat_real = 0 ! Number of real properties per material mixed_int = 0 ! Number of integer properties per mixed_bc mixed_real = 0 ! Number of real properties per mixed_bc pt_int = 0 ! Number of integer properties per node pt_real = 0 ! Number of real properties per node reals = 0 ! Number of miscellaneous real properties seg_int = 0 ! Number of integer properties per segment seg_real = 0 ! Number of real properties per segment ! Logicals ANISOTROPIC = .FALSE. ! Need anisotropic constitutive law AVERAGE_MASS = .FALSE. ! Average consistent and diagonal mass AXISYMMETRIC = .FALSE. ! Axisymmetrix solid, x=radius BAR_CHART = .FALSE. ! No bar charts in output BUCKLING = .FALSE. ! No buckling eigenanalysis CONVECTION = .FALSE. ! No constant convection mixed bc CONVECT_VARY = .FALSE. ! Convection, different on each mixed seg CONSTANT_J = .FALSE. ! Element Jacobian is not constant DATA_CHECK_ONLY = .FALSE. ! Stop after reading the input data DEBUG_ADAPT = .FALSE. ! No debug print of adaptive items DEBUG_ASSEMBLY = .FALSE. ! No debug print of assembled SS, CC DEBUG_B = .FALSE. ! No debug print of B matrix DEBUG_E = .FALSE. ! No debug print of E matrix DEBUG_EL_COL = .FALSE. ! No debug print of element column matrix DEBUG_EL_SQ = .FALSE. ! No debug print of element square matrix DEBUG_EL_POST = .FALSE. ! No debug print of saving post-processing info DEBUG_EL_TYPE = .FALSE. ! No debug print of element type properties DEBUG_ERROR = .FALSE. ! No debug print of error estimator DEBUG_EXACT = .FALSE. ! No debug print of exact case items DEBUG_INERTIA = .FALSE. ! No debug print of geometric inertia DEBUG_INCLUDE = .FALSE. ! No debug print of include files DEBUG_MIX_SQ = .FALSE. ! No debug print of mixed BC matrices DEBUG_POST_EL = .FALSE. ! No debug print of post-processing an element DEBUG_PROPERTY = .FALSE. ! No debug print of any GET property DEBUG_SCP = .FALSE. ! No debug print of SCP averaging DEBUG_SEG_COL = .FALSE. ! No debug print of segment column matrix DEBUG_UNITS = .FALSE. ! No debug print of file unit numbers DIAGONAL_MASS = .FALSE. ! Diagonalize mass matrix DOF_VARY = .FALSE. ! Number of nodal dof can vary DYNAMIC = .FALSE. ! Problem has 2-nd derivative in time EIGEN = .FALSE. ! Problem is eigen analysis EIGEN_ROOT = .FALSE. ! Report square root of eigenvalue EL_DAMPING = .FALSE. ! Element damping sq matrix is defined FLUX_NORMAL = .FALSE. ! Problem has constant normal flux FORCE_BODY = .FALSE. ! No body force per unit volume input FOURIER_LAW = .FALSE. ! Fourier flux law, q = -E*Del_U FREQ_CONTROL = .FALSE. ! No force Sin or Cos frequency given FROM_REST = .FALSE. ! Problem starts from rest, v=a=0 F90 = .FALSE. ! F90 compiler was used (for maxloc) F95 = .FALSE. ! F95 compiler was used (for maxloc) GEOMETRIC_K = .FALSE. ! Form geometric stiffness matrix GET_DELETED_EL = .FALSE. ! Read list of deleted elements GET_2ND_DERIV = .FALSE. ! Build 2nd derivatives of H array GLOBAL_PROPERTY = .FALSE. ! Was a global property set by keyword GRAD_BASE_ERROR = .FALSE. ! Base error estimates on gradients only HAS_BODY_FORCE = .FALSE. ! No body force per unit volume input HEAVYSIDE = .FALSE. ! No Heavyside step up or down input HUGHES_METHOD = .FALSE. ! Use Hughes dynamic time integration INCREMENTAL = .FALSE. ! Use an incremental solution INITIAL_FUNCTION = .FALSE. ! Initial condition is spatial function INPUT_ONLY = .FALSE. ! Just read and validate data IS_ELEMENT = .FALSE. ! Element class unknown IS_FLUX_SEG = .FALSE. ! Element class unknown IS_MIXED_BC = .FALSE. ! Element class unknown IS_NULL = .FALSE. ! Solution values will be zero JACOBI = .FALSE. ! Full matrix general eigenvalue solution LIST_EL_TO_EL = .FALSE. ! List element neighbors LIST_QP_FLUX = .FALSE. ! Omit list of flux at every qp in element NEWMARK_METHOD = .TRUE. ! Use Newmark time integration method NODAL_TAU = .FALSE. ! Use nodal basis SUPG method NO_ERROR_EST = .FALSE. ! Expect to do SCP error estimates NO_PRINTING = .FALSE. ! Do not print results (except *.tmp) NO_PROPERTIES = .FALSE. ! Expect property input NO_SCP_AVE = .FALSE. ! Do not save data for SCP averages PLANE_STRAIN = .FALSE. ! Use 2-D plane stress assumption POLAR_COORD = .FALSE. ! Four dimensional polar space POST_EL = .FALSE. ! No user standard element post-process POST_MIXED = .FALSE. ! No user mixed element post-process PROPERTY_FIRST = .FALSE. ! Read all properties after remarks PT_EL_LIST = .FALSE. ! List all elements at each node PT_EL_SUMS = .FALSE. ! List element count at each node PT_MASSES = .FALSE. ! Point masses will be input PT_DAMPING = .FALSE. ! Point dampners will be input RAMP_CONTROL = .FALSE. ! No load ramp up or down given SAVE_NEW_MESH = .FALSE. ! Save new sizes for adaptive mesher SAVE_1248 = .FALSE. ! Save after steps 1, 2, 4, 8, 16, SCALE_FREQ = .FALSE. ! Time scale based on frequency SCALE_LATE = .FALSE. ! Time scale starts after 0 SCALE_LOADS = .FALSE. ! Scale loads in iteration or time SCALE_RAMP = .FALSE. ! Time scale is a ramp to/from 1/0 SCALE_STEP = .FALSE. ! Time scale is a step to/from 1/0 SCALE_ZERO = .FALSE. ! Time scale starts at 0 SCP_CENTER_ONLY = .FALSE. ! Average center node or element only SCP_DOF = .FALSE. ! Average patch degrees of freedom only SCP_E_AT_NP = .FALSE. ! Variable E defined at nodal points !b SCP_E_AT_QP = .FALSE. ! Variable E defined at quadrature pts !b SCP_E_SMOOTH = .FALSE. ! E smoothly varies over element SCP_2ND_DERIV = .FALSE. ! Do not fit second derivatives in SCP SCP_3RD_DERIV = .FALSE. ! Do not fit third derivatives in SCP SCP_NEIGH_EL = .TRUE. ! Patch includes all elements SCP_NEIGH_FACE = .FALSE. ! Patch includes facing elements SCP_NEIGH_PT = .FALSE. ! Patch includes elements adjacent to node SCP_ONLY_ONCE = .TRUE. ! Use only averages from 1st el in patch SCP_REC_NUMB_ALLOC = .FALSE. ! SCP memory allocation SCP_RECORDS_SAVED = .FALSE. ! SCP record status SKIP_ERROR = .FALSE. ! Skip the SCP error estimates SPACE_TIME = .FALSE. ! A unstructured space-time formulation S_T_CONTINUOUS = .FALSE. ! A space-time_slab continuous formulation STATIC = .TRUE. ! Time independent problem SUPG = .FALSE. ! Use the SUPG upwind method SYMMETRIC = .TRUE. ! Symmetric skyline equations TAU_BOX = .FALSE. ! Akin bounding box method for SUPG Tau TAU_GEOM = .FALSE. ! Akin geometry method for SUPG Tau TAU_PEC = .FALSE. ! Peclet 1-D opt method for SUPG Tau TAU_QP = .FALSE. ! Form Tau at quadrature points TAU_SHOW = .FALSE. ! Show all Tau values in first element TAU_S1 = .FALSE. ! Tezduyar S_1 method for SUPG Tau TAU_RGN = .FALSE. ! Tezduyar RGN method for SUPG Tau TAU_UGN = .FALSE. ! Tezduyar UGN method for SUPG Tau TAU_VOL = .FALSE. ! Akin volume method for SUPG Tau TIME_SLAB = .FALSE. ! A space-time_slab formulation TOUCH = .FALSE. ! Dummy to suppress compile warnings TRANSIENT = .FALSE. ! Problem has 1-st derivative in time TYPES_ALLOCATED = .FALSE. ! Element types control data allocated UPWIND = .FALSE. ! Upwind the advection term USE_EXACT = .FALSE. ! Use exact solution in error est, i/o USE_EXACT_BC = .FALSE. ! Use exact solution in essential bc USE_EXACT_FLUX = .FALSE. ! Use exact fluxes in error est, i/o USE_EXACT_ROBIN = .FALSE. ! Use exact solution Robin BC terms USE_EXACT_SOURCE = .FALSE. ! Use exact solution source term USER_LOGIC = .FALSE. ! Logical variable for user control VARY_E = .FALSE. ! Vary the E matrix inside each element WILSON_METHOD = .FALSE. ! Use Wilson time integration method END SUBROUTINE SET_CONTROL_DEFAULTS SUBROUTINE ALLOCATE_TYPE_CONTROLS IF ( .NOT. TYPES_ALLOCATED ) THEN IF ( N_L_TYPE == 1 ) THEN PRINT *, 'WARNING: ELEMENT TYPES ACTIVE FOR ONE TYPE' PRINT *, 'SET KEYWORD el_types > 1' N_WARN = N_WARN + 1 END IF ALLOCATE ( TYPE_NODES (N_L_TYPE) ) ALLOCATE ( TYPE_GAUSS (N_L_TYPE) ) ALLOCATE ( TYPE_GEOM (N_L_TYPE) ) ALLOCATE ( TYPE_PARM (N_L_TYPE) ) ALLOCATE ( TYPE_SHAPE (N_L_TYPE) ) ALLOCATE ( TYPE_FACES (N_L_TYPE) ) ALLOCATE ( TYPE_SONS (N_L_TYPE) ) ALLOCATE ( TYPE_USER (N_L_TYPE) ) TYPES_ALLOCATED = .TRUE. TYPE_NODES = 0 ; TYPE_GAUSS = 0 ; TYPE_GEOM = 0 ; TYPE_PARM = 0 TYPE_SHAPE = 0 ; TYPE_FACES = 0 ; TYPE_SONS = 0 ; TYPE_USER = 0 ELSE PRINT *,'WARNING: Element types controls already allocated' N_WARN = N_WARN + 1 END IF END SUBROUTINE ALLOCATE_TYPE_CONTROLS SUBROUTINE ALLOCATE_OUTPUT_NAMES ! ASSIGN NAMES TO VARIOUS OUTPUT LISTS CALL ALLOCATE_DOF_NAMES ! and INPUT_NAMES CALL ALLOCATE_XYZ_NAMES ! and NO_NAME CALL ALLOCATE_FLUX_NAMES CALL ALLOCATE_EXACT_NAMES CALL ALLOCATE_SCP_AVE_NAMES IF ( N_2_DER > 0 ) CALL ALLOCATE_DERIV_2_NAMES IF ( N_3_DER > 0 ) CALL ALLOCATE_DERIV_3_NAMES END SUBROUTINE ALLOCATE_OUTPUT_NAMES SUBROUTINE ALLOCATE_SCP_AVE_NAMES ! ASSIGN NAMES TO ITEMS AVERAGED WITH SCP PROCESS CHARACTER(len=1) :: J_CHAR ; INTEGER :: J ! Loop ALLOCATE ( AVE_NAME (SCP_FIT), AVE_SHORT (SCP_FIT) ) IF ( SCP_FIT > 9 ) THEN PRINT *,'WARNING, ALLOCATE_AVE_NAMES: MORE THAN 9 GENERAL DOF' N_WARN = N_WARN + 1 END IF ! Invalid name creation DO J = 1, SCP_FIT WRITE (J_CHAR, '(I1)') J AVE_NAME (J) = ' SCP_AVE_' // J_CHAR // ',' AVE_SHORT (J) = ' AVE_' // J_CHAR // ',' END DO END SUBROUTINE ALLOCATE_SCP_AVE_NAMES SUBROUTINE ALLOCATE_DOF_NAMES ! ASSIGN NAMES TO DEGREES OF FREEDOM CHARACTER(len=1) :: J_CHAR ; INTEGER :: J ! Loop ALLOCATE ( DOF_NAME (N_G_DOF), DOF_SHORT (N_G_DOF) ) ALLOCATE ( INPUT_NAME (N_G_DOF) ) IF ( N_G_DOF > 9 ) THEN PRINT *,'WARNING, ALLOCATE_DOF_NAMES: MORE THAN 9 GENERAL DOF' N_WARN = N_WARN + 1 END IF ! Invalid name creation DO J = 1, N_G_DOF WRITE (J_CHAR, '(I1)') J DOF_NAME (J) = ' DOF_' // J_CHAR // ',' DOF_SHORT (J) = ' D_' // J_CHAR // ',' INPUT_NAME (J) = ' IN_' // J_CHAR // ',' END DO END SUBROUTINE ALLOCATE_DOF_NAMES SUBROUTINE ALLOCATE_EXACT_NAMES ! ASSIGN NAMES TO EXACT DEGREES OF FREEDOM CHARACTER(len=1) :: J_CHAR ; INTEGER :: J, LARGE ! Loop LARGE = MAX (N_G_DOF, N_R_B, N_2_DER, N_3_DER, SCP_FIT) !b ALLOCATE ( EXACT_NAME (N_G_DOF) ) ALLOCATE ( EXACT_NAME (LARGE), EXACT_SHORT (LARGE) ) !b IF ( LARGE > 9 ) THEN PRINT *,'WARNING, ALLOCATE_EXACT_NAMES: MORE THAN 9 ENTRIES' N_WARN = N_WARN + 1 END IF ! Invalid name creation !b DO J = 1, N_G_DOF DO J = 1, LARGE !b WRITE (J_CHAR, '(I1)') J EXACT_NAME (J) = ' EXACT' // J_CHAR // ',' EXACT_SHORT (J) = ' EX_' // J_CHAR // ',' END DO END SUBROUTINE ALLOCATE_EXACT_NAMES SUBROUTINE ALLOCATE_FLUX_NAMES ! ASSIGN NAMES TO COMPUTED GENERALIZED FLUX COMPONENTS ! (NOT INPUT SURFACE FLUX COMPONENTS) CHARACTER(len=1) :: J_CHAR CHARACTER(len=2) :: J_CHAR2 INTEGER :: J, HIGHER, LARGE ! Loop HIGHER = N_R_B * N_SPACE LARGE = MAX ( HIGHER, (N_R_B + 2)) + N_2_DER ALLOCATE ( FLUX_NAME (LARGE), FLUX_SHORT (LARGE) ) !b IF ( LARGE > 9 ) THEN !b PRINT *,'WARNING, ALLOCATE_FLUX_NAMES: NEED MORE CHARACTERS' !b PRINT *,'N_R_B, N_SPACE, N_2_DER ', N_R_B, N_SPACE, N_2_DER !b N_WARN = N_WARN + 1 !b END IF ! Invalid name creation DO J = 1, LARGE IF ( J <= 9 ) THEN WRITE (J_CHAR, '(I1)') J FLUX_NAME (J) = ' FLUX_' // J_CHAR // ',' FLUX_SHORT (J) = ' FX_' // J_CHAR // ',' ELSE WRITE (J_CHAR2, '(I2)') J FLUX_NAME (J) = ' FLUX_' // J_CHAR2 // ',' FLUX_SHORT (J) = ' FX_' // J_CHAR2 // ',' END IF END DO END SUBROUTINE ALLOCATE_FLUX_NAMES SUBROUTINE ALLOCATE_DERIV_2_NAMES ! ASSIGN NAMES TO COMPUTED GENERALIZED DERIV_2 COMPONENTS ! (NOT INPUT SURFACE DERIV_2 COMPONENTS) CHARACTER(len=1) :: J_CHAR ; INTEGER :: J, HIGHER, LARGE ! Loop HIGHER = N_R_B * N_SPACE LARGE = MAX ( HIGHER, (N_2_DER + 2) ) ALLOCATE ( DERIV_2_NAME (LARGE), DERIV_2_SHORT (LARGE) ) IF ( LARGE > 9 ) THEN PRINT *,'WARNING, ALLOCATE_DERIV_2_NAMES: MORE THAN 9 COMPONENTS' N_WARN = N_WARN + 1 END IF ! Invalid name creation DO J = 1, LARGE WRITE (J_CHAR, '(I1)') J DERIV_2_NAME (J) = ' DERIV_2_' // J_CHAR // ',' DERIV_2_SHORT (J) = ' DRV_2_' // J_CHAR // ',' END DO IF ( N_G_DOF == 1 .AND. N_SPACE == 2 ) THEN ! Special names !'123456789012' DERIV_2_NAME (1) = ' U_xx '; DERIV_2_SHORT (1) = ' U_xx ' DERIV_2_NAME (2) = ' U_yx '; DERIV_2_SHORT (2) = ' U_yx ' DERIV_2_NAME (3) = ' U_xy '; DERIV_2_SHORT (3) = ' U_xy ' DERIV_2_NAME (4) = ' U_yy '; DERIV_2_SHORT (4) = ' U_yy ' END IF END SUBROUTINE ALLOCATE_DERIV_2_NAMES SUBROUTINE ALLOCATE_DERIV_3_NAMES ! ASSIGN NAMES TO COMPUTED GENERALIZED DERIV_2 COMPONENTS ! (NOT INPUT SURFACE DERIV_2 COMPONENTS) CHARACTER(len=1) :: J_CHAR ; INTEGER :: J ! Loop ALLOCATE ( DERIV_3_NAME (N_3_DER + 2) ) !bIF ( (N_3_DER + 2) > 9 ) THEN !bPRINT *,'WARNING, ALLOCATE_DERIV_3_NAMES: > 9 COMPONENTS' !bN_WARN = N_WARN + 1 !bEND IF ! Invalid name creation DO J = 1, (N_3_DER + 2) WRITE (J_CHAR, '(I1)') J DERIV_3_NAME (J) = ' DERIV_3_' // J_CHAR // ',' END DO END SUBROUTINE ALLOCATE_DERIV_3_NAMES SUBROUTINE ALLOCATE_XYZ_NAMES ! ASSIGN NAMES TO SPACIAL COORDINATES CHARACTER(len=12), PARAMETER :: BLANK = ' ' ! blank name CHARACTER(len=1) :: J_CHAR ; INTEGER :: J ! Loop ALLOCATE ( XYZ_NAME (N_SPACE), NO_NAME (N_SPACE) ) ALLOCATE ( XYZ_SHORT (N_SPACE), NO_SHORT (N_SPACE) ) NO_NAME = BLANK ; NO_SHORT = BLANK SELECT CASE (N_SPACE) CASE (1:3) ; XYZ_NAME (1:N_SPACE) = RECTANGULAR (1:N_SPACE) XYZ_SHORT (1:N_SPACE) = RECT_SHORT (1:N_SPACE) IF (AXISYMMETRIC) THEN XYZ_NAME (1:N_SPACE) = CYLINDRICAL (1:N_SPACE) XYZ_SHORT (1:N_SPACE) = CYL_SHORT (1:N_SPACE) END IF ! axisymmetric CASE ( 4 ) ; XYZ_NAME (1:N_SPACE) = POLAR_SPACE (1:N_SPACE) XYZ_SHORT (1:N_SPACE) = POLAR_SHORT (1:N_SPACE) CASE DEFAULT ! Build the names IF ( N_SPACE > 9 ) THEN PRINT *,'WARNING, ALLOCATE_XYZ_NAMES: MORE THAN 9 SPACES ACTIVE' N_WARN = N_WARN + 1 END IF ! Invalid name creation DO J = 1, N_SPACE WRITE (J_CHAR, '(I1)') J XYZ_NAME (J) = ' X_' // J_CHAR // ' COORD, ' XYZ_SHORT (J) = ' X_' // J_CHAR // ', ' END DO END SELECT END SUBROUTINE ALLOCATE_XYZ_NAMES SUBROUTINE INITIALIZE_SCP_RECORDS () IF ( U_SCPR > 0 .AND. N_QP > 0 .AND. (.NOT. SCP_REC_NUMB_ALLOC) ) THEN ALLOCATE ( SCP_RECORD_NUMBER (N_ELEMS, N_QP) ) SCP_REC_NUMB_ALLOC = .TRUE. SCP_RECORDS_SAVED = .FALSE. SCP_RECORD_NUMBER = 0 SCP_PAD = EST_SCP_PAD_SIZE () SCP_RECL = SET_SCP_RECORD_SIZE () ! OPEN DIRECT (RANDOM) ACCESS SCP FILE IF ( SCP_RECL > 0 ) THEN ! VALID RECORD LENGTH OPEN (U_SCPR, STATUS = 'SCRATCH', FORM = 'UNFORMATTED', & ACCESS = 'DIRECT', RECL = SCP_RECL, ACTION ='READWRITE') !b ACCESS = 'DIRECT', RECL = SCP_RECL, ACTION ='READWRITE',& !b FILE = 'U_SCPR_TMP.BIN') PRINT *, ' ' ; PRINT *, 'SUPER_CONVERGENT PATCH FILE OPENED' ELSE ! ISSUE WARNING PRINT *, 'WARNING, SCPR RECORD LENGTH NOT SET' PRINT *, 'USE SET_SCPR_RECORD_LENGTH' U_SCPR = 0 ! TURN OFF SCP OPTION N_WARN = N_WARN + 1 ! INCREMENT WARNING END IF ! VALID DIRECT RECORD LENGTH PRINT *, 'SUPER_CONVERGENT PATCH RECORDS INITIALIZED' IF ( SCP_RECL < 1 ) THEN ! INVALID RECORD LENGTH PRINT *, 'WARNING, SCPR RECORD LENGTH NOT SET' PRINT *, 'USE SET_SCPR_RECORD_LENGTH IN ELEM_SQ_MATRIX' N_WARN = N_WARN + 1 ! INCREMENT WARNING END IF ! SCP CHECK ELSE IF ( THIS_ITER == 1 ) THEN !b PRINT *, 'WARNING: NO SUPER_CONVERGENT PATCH RECORDS AVAILABLE' N_WARN = N_WARN + 1 ! INCREMENT WARNING END IF ! iterations possible END IF ! RECORDS NEEDED END SUBROUTINE INITIALIZE_SCP_RECORDS SUBROUTINE DEALLOCATE_SCP_RECORDS DEALLOCATE ( SCP_RECORD_NUMBER ) SCP_REC_NUMB_ALLOC = .FALSE. SCP_RECORDS_SAVED = .FALSE. END SUBROUTINE DEALLOCATE_SCP_RECORDS FUNCTION SET_SCP_RECORD_SIZE () RESULT (VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * * ! SET RECL TO STORE: XYZ, E, STRESS AT EACH QP ! * * * * * * * * * * * * * * * * * * * * * * * * IMPLICIT NONE INTEGER :: VALUE REAL(DP) :: FAKE (N_SPACE + N_R_B*N_R_B + SCP_FIT + SCP_PAD) INQUIRE ( IOLENGTH = VALUE ) FAKE IF ( I_BUG > 0 ) PRINT *, 'SCP RECORD SIZE SET TO = ', VALUE END FUNCTION SET_SCP_RECORD_SIZE SUBROUTINE UPDATE_SCP_STATUS ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * ! DID LIST_ELEM_GRADIENTS OR LIST_ELEM_FLUXES SAVE DATA ! FROM A CALL IN POST_PROCESS_ELEM ? ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * IF ( RECORD_NUMBER > 0 ) THEN ! DATA WERE SAVED SCP_RECORDS_SAVED = .TRUE. ELSE SCP_RECORDS_SAVED = .FALSE. END IF END SUBROUTINE UPDATE_SCP_STATUS FUNCTION EST_SCP_PAD_SIZE () RESULT (VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * * ! ESTIMATE PAD SIZE TO USE FULL SECTORS FOR ! SCP RECL TO STORE: XYZ, E, STRESS ! * * * * * * * * * * * * * * * * * * * * * * * * IMPLICIT NONE INTEGER :: LENGTH, VALUE REAL(DP) :: FAKE (N_SPACE + N_R_B*N_R_B + SCP_FIT) INTEGER, PARAMETER :: SECTOR_SIZE = 512 INQUIRE ( IOLENGTH = LENGTH ) FAKE ! HARDWARE RECORD LENGTH VALUE = LENGTH IF ( LENGTH < SECTOR_SIZE ) VALUE = SECTOR_SIZE - LENGTH IF ( I_BUG > 0 ) PRINT *, 'RECOMMENDED SCP PAD SIZE = ', VALUE VALUE = 0 !b force to zero in development END FUNCTION EST_SCP_PAD_SIZE SUBROUTINE TURN_OFF_PLOTS ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * ! TURN OFF SAVE OF FORMATTED RESULTS (AND MATLAB PLOTS) ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * U_PLT1 = 0 ; U_PLT2 = 0 ; U_PLT3 = 0; U_PLT4 = 0 ; U_PLT5 = 0 ; U_PLT6 = 0 U2_PLT1 = 0 ; U2_PLT2 = 0 ; U2_PLT3 = 0; U2_PLT4 = 0 END SUBROUTINE TURN_OFF_PLOTS SUBROUTINE TURN_ON_BAR_CHART_DIST NO_DIST = 0 ! SCALE DISTANCE BETWEEN NODES IN CHART END SUBROUTINE TURN_ON_BAR_CHART_DIST subroutine set_N_FILE3 (k) integer, intent(in) :: k N_FILE3 = k if ( N_FILE3 > 0 ) open (N_FILE3, file="N_FILE3.FIL", & status="unknown") end subroutine set_N_FILE3 SUBROUTINE SET_AXISYMMETRIC_TO_TRUE AXISYMMETRIC = .TRUE. END SUBROUTINE SET_AXISYMMETRIC_TO_TRUE !b SUBROUTINE ALLOCATE_SYSTEM_ARRAYS !b ALLOCATE ( BODY_F (N_SPACE) ) ! will be body !b ALLOCATE ( I_BC (MAX_NP) ) !b ALLOCATE ( NODES (L_S_TOT, NOD_PER_EL) ) !b ALLOCATE ( N_RES_EQ (MAX_TYP) ) !b ALLOCATE ( X (MAX_NP, N_SPACE) ) !b if (.not. ALLOCATED (NODES )) stop 'NODES not allocated' !b if (.not. ALLOCATED (N_RES_EQ)) stop 'N_RES_EQ not allocated' !b if (.not. ALLOCATED (I_BC )) stop 'I_BC not allocated' !b if (.not. ALLOCATED (X )) stop 'X not allocated' !b END SUBROUTINE ALLOCATE_SYSTEM_ARRAYS FUNCTION GET_THIS_ELEMENT_NUMBER () RESULT (IE) INTEGER :: IE IF ( THIS_EL < 1 ) THEN PRINT *,'ERROR: GET_THIS_ELEMENT_NUMBER THIS_EL NOT SET' PRINT *,'SET_THIS_ELEMENT_NUMBER (IE) HAS NOT BEEN USED' STOP 'ERROR: GET_THIS_ELEMENT_NUMBER THIS_EL NOT SET' END IF IE = THIS_EL END FUNCTION GET_THIS_ELEMENT_NUMBER FUNCTION GET_THIS_MIXED_SEG_NUMBER () RESULT (IE) INTEGER :: IE IF ( IS_MIX < 1 ) THEN PRINT *,'ERROR: GET_THIS_MIXED_SEG_NUMBER IS_MIX NOT SET' PRINT *,'SET_THIS_MIXED_SEG_NUMBER (IE) HAS NOT BEEN USED' STOP 'ERROR: GET_THIS_MIXED_SEG_NUMBER IS_MIX NOT SET' END IF IE = IS_MIX END FUNCTION GET_THIS_MIXED_SEG_NUMBER FUNCTION GET_THIS_FLUX_SEG_NUMBER () RESULT (IE) INTEGER :: IE IF ( IS_SEG < 1 ) THEN PRINT *,'ERROR: GET_THIS_FLUX_SEG_NUMBER IS_SEG NOT SET' PRINT *,'SET_THIS_FLUX_SEG_NUMBER (IE) HAS NOT BEEN USED' STOP 'ERROR: GET_THIS_FLUX_SEG_NUMBER IS_SEG NOT SET' END IF IE = IS_SEG END FUNCTION GET_THIS_FLUX_SEG_NUMBER FUNCTION GET_THIS_ITERATION_NUMBER () RESULT (ITER) INTEGER :: ITER ITER = THIS_ITER END FUNCTION GET_THIS_ITERATION_NUMBER SUBROUTINE SET_THIS_ELEMENT_NUMBER (IE) INTEGER, INTENT (IN) :: IE THIS_EL = IE END SUBROUTINE SET_THIS_ELEMENT_NUMBER SUBROUTINE CITE_INPUT_LINE ! LIST THE CURRENT INPUT LINE NUMBER PRINT *, 'NOTE: CURRENTLY AT DATA INPUT LINE NUMBER ', & LINE_IN END SUBROUTINE CITE_INPUT_LINE SUBROUTINE INCREMENT_INPUT_LINE ! INCREMENT THE CURRENT INPUT LINE NUMBER LINE_IN = LINE_IN + 1 !b print *, 'LINE_IN ', LINE_IN !b END SUBROUTINE INCREMENT_INPUT_LINE End Module System_Constants Module Elem_Type_Data Use System_Constants IMPLICIT NONE INTEGER :: LAST_LT, LT_FREE INTEGER :: LT_N, LT_QP, LT_GEOM, LT_PARM INTEGER :: LT_FACES, LT_SONS, LT_SHAP, LT_USER LOGICAL :: TYPE_APLY_ALLOC = .FALSE. LOGICAL :: TYPE_DATA_ALLOC = .FALSE. LOGICAL :: TYPE_DOF_ALLOC = .FALSE. LOGICAL :: TYPE_EQS_ALLOC = .FALSE. LOGICAL :: TYPE_GAUS_ALLOC = .FALSE. LOGICAL :: TYPE_GEOM_ALLOC = .FALSE. LOGICAL :: TYPE_NTRP_ALLOC = .FALSE. LOGICAL :: TYPE_SCAL_ALLOC = .FALSE. LOGICAL :: TYPE_TOPO_ALLOC = .FALSE. INTEGER, ALLOCATABLE :: L_TYPE (:) ! type number of each elem INTEGER, ALLOCATABLE :: LT_DATA (:,:) ! constants defining a type ! type quadrature points REAL(DP), ALLOCATABLE :: WT (:), PT (:, :) ! quadratures ! DOF's and derivatives at type quadrature points REAL(DP), ALLOCATABLE :: V_QP (:, :), DLV_QP (:, :, :) ! vec filled REAL(DP), ALLOCATABLE :: V (:), DLV (:, :) ! vec at pt ! geometry and derivatives at type quadrature points REAL(DP), ALLOCATABLE :: G_QP (:, :), DLG_QP (:, :, :) ! geom filled REAL(DP), ALLOCATABLE :: G (:), DLG (:, :) ! geom at pt ! scalar and derivatives at type quadrature points REAL(DP), ALLOCATABLE :: H_QP (:, :), DLH_QP (:, :, :) ! scalar filled REAL(DP), ALLOCATABLE :: H (:), DLH (:, :) ! scalar at pt REAL(DP), ALLOCATABLE :: D2LH (:, :) ! scalar at pt ! Basic type geometry arrays INTEGER, ALLOCATABLE :: ELEM_NODES (:) ! node numbers REAL(DP), ALLOCATABLE :: COORD (:, :) ! all spatial coordinates REAL(DP), ALLOCATABLE :: GEOMETRY (:, :) ! geometric coordinates ! Application specific element type arrays INTEGER, ALLOCATABLE :: INDEX (:) ! equation numbers REAL(DP), ALLOCATABLE :: D (:) ! element solution vector REAL(DP), ALLOCATABLE :: C (:) ! element source vector REAL(DP), ALLOCATABLE :: S (:, :) ! element square matrix REAL(DP), ALLOCATABLE :: DAMP (:, :) ! element damping matrix REAL(DP), ALLOCATABLE :: EL_M (:, :) ! element mass matrix REAL(DP), ALLOCATABLE :: DIAG_M (:) ! diagonal mass matrix contains SUBROUTINE SET_ELEM_TYPE_DEFAULTS ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * ! Set element type defaults to override with keyword inputs ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * ! Integers LAST_LT = 0 ! Previous element type number LT_FACES = 0 ! Number of element type faces LT_FREE = 0 ! Number of degrees-of-freedom for element type LT_GEOM = 0 ! Number of element type geometric nodes LT_N = 0 ! Number of nodes per element type LT_PARM = 0 ! Dimension of element type parametric space LT_QP = 0 ! Number of quadrature points for element type LT_SHAP = 0 ! Element type shape flag number LT_SONS = 0 ! Number of element type sons, if divided LT_USER = 0 ! Application dependent optional user item ! Logicals TYPE_APLY_ALLOC = .FALSE. ! Overall status of element type TYPE_DATA_ALLOC = .FALSE. ! Overall status of master system data TYPE_DOF_ALLOC = .FALSE. ! Vector interpolation & gradients TYPE_EQS_ALLOC = .FALSE. ! Element matrices and assembly index TYPE_GAUS_ALLOC = .FALSE. ! Quadrature arrays TYPE_GEOM_ALLOC = .FALSE. ! Geometry interpolation & gradients TYPE_NTRP_ALLOC = .FALSE. ! Overall status of interpolations TYPE_SCAL_ALLOC = .FALSE. ! Scalar interpolation & gradients TYPE_TOPO_ALLOC = .FALSE. ! Node (topology) arrays END SUBROUTINE SET_ELEM_TYPE_DEFAULTS SUBROUTINE LIST_TYPE_ALLOC_STATUS PRINT *, ' ' PRINT *, '** TYPE ALLOCATION STATUS **' PRINT *, 'TYPE_APLY_ALLOC = ', TYPE_APLY_ALLOC PRINT *, 'TYPE_DATA_ALLOC = ', TYPE_DATA_ALLOC PRINT *, 'TYPE_DOF_ALLOC = ', TYPE_DOF_ALLOC PRINT *, 'TYPE_EQS_ALLOC = ', TYPE_EQS_ALLOC PRINT *, 'TYPE_GAUS_ALLOC = ', TYPE_GAUS_ALLOC PRINT *, 'TYPE_GEOM_ALLOC = ', TYPE_GEOM_ALLOC PRINT *, 'TYPE_NTRP_ALLOC = ', TYPE_NTRP_ALLOC PRINT *, 'TYPE_SCAL_ALLOC = ', TYPE_SCAL_ALLOC PRINT *, 'TYPE_TOPO_ALLOC = ', TYPE_TOPO_ALLOC END SUBROUTINE LIST_TYPE_ALLOC_STATUS SUBROUTINE ALLOCATE_TYPE_APPLICATION IF ( .NOT. TYPE_TOPO_ALLOC ) THEN ALLOCATE ( ELEM_NODES (LT_N) ) ! TOPOLOGY ALLOCATE ( COORD (LT_N, N_SPACE) ) ALLOCATE ( GEOMETRY (LT_GEOM, N_SPACE) ) TYPE_TOPO_ALLOC = .TRUE. IF ( I_BUG > 0 ) PRINT *, 'ALLOCATED TYPE TOPOLOGY' ELSE PRINT *, 'WARNING: TYPE TOPO ALREADY ALLOCATED' N_WARN = N_WARN + 1 ! INCREMENT WARNING END IF IF ( .NOT. TYPE_EQS_ALLOC ) THEN ALLOCATE ( INDEX (LT_FREE) ) ! ELEMENT EQUATIONS ALLOCATE ( D (LT_FREE) ) ALLOCATE ( C (LT_FREE) ) ALLOCATE ( S (LT_FREE, LT_FREE) ) ALLOCATE ( DAMP (LT_FREE, LT_FREE) ) ALLOCATE ( EL_M (LT_FREE, LT_FREE) ) ALLOCATE ( DIAG_M (LT_FREE) ) TYPE_EQS_ALLOC = .TRUE. IF ( I_BUG > 0 ) PRINT *, 'ALLOCATED TYPE EQUATIONS' ELSE PRINT *, 'WARNING: TYPE EQS ALREADY ALLOCATED' N_WARN = N_WARN + 1 ! INCREMENT WARNING END IF TYPE_APLY_ALLOC = .TRUE. END SUBROUTINE ALLOCATE_TYPE_APPLICATION SUBROUTINE DEALLOCATE_TYPE_APPLICATION IF ( TYPE_EQS_ALLOC ) THEN DEALLOCATE ( DIAG_M ) DEALLOCATE ( EL_M ) DEALLOCATE ( DAMP ) DEALLOCATE ( S ) DEALLOCATE ( C ) DEALLOCATE ( D ) DEALLOCATE ( INDEX ) TYPE_EQS_ALLOC = .FALSE. IF ( I_BUG > 0 ) PRINT *, 'DEALLOCATED TYPE EQUATIONS' ELSE PRINT *, 'WARNING: TYPE EQS ALREADY DEALLOCATED' N_WARN = N_WARN + 1 ! INCREMENT WARNING END IF IF ( TYPE_APLY_ALLOC ) THEN DEALLOCATE ( GEOMETRY ) DEALLOCATE ( COORD ) DEALLOCATE ( ELEM_NODES ) TYPE_TOPO_ALLOC = .FALSE. IF ( I_BUG > 0 ) PRINT *, 'DEALLOCATED TYPE TOPOLOGY' ELSE PRINT *, 'WARNING: TYPE APP ALREADY DEALLOCATED' N_WARN = N_WARN + 1 ! INCREMENT WARNING END IF TYPE_APLY_ALLOC = .FALSE. END SUBROUTINE DEALLOCATE_TYPE_APPLICATION SUBROUTINE ALLOCATE_TYPE_QUADRATURE IF ( .NOT. TYPE_GAUS_ALLOC ) THEN ALLOCATE ( WT (LT_QP) ) ALLOCATE ( PT (LT_PARM, LT_QP) ) TYPE_GAUS_ALLOC = .TRUE. IF ( I_BUG > 0 ) PRINT *, 'ALLOCATED TYPE QUADRATURE' ELSE PRINT *, 'WARNING: TYPE QP ALREADY ALLOCATED' N_WARN = N_WARN + 1 ! INCREMENT WARNING END IF END SUBROUTINE ALLOCATE_TYPE_QUADRATURE SUBROUTINE DEALLOCATE_TYPE_QUADRATURE IF ( TYPE_GAUS_ALLOC ) THEN DEALLOCATE ( PT ) DEALLOCATE ( WT ) TYPE_GAUS_ALLOC = .FALSE. IF ( I_BUG > 0 ) PRINT *, 'DEALLOCATED TYPE QUADRATURE' ELSE PRINT *, 'WARNING: TYPE QP ALREADY DEALLOCATED' N_WARN = N_WARN + 1 ! INCREMENT WARNING END IF END SUBROUTINE DEALLOCATE_TYPE_QUADRATURE SUBROUTINE ALLOCATE_TYPE_INTERPOLATIONS IF ( .NOT. TYPE_GEOM_ALLOC ) THEN CALL ALLOCATE_TYPE_QUADRATURE ALLOCATE ( G (LT_GEOM) ) ! GEOMETRY ALLOCATE ( G_QP (LT_GEOM, LT_QP) ) ALLOCATE ( DLG (LT_PARM, LT_GEOM) ) ALLOCATE ( DLG_QP (LT_PARM, LT_GEOM, LT_QP) ) TYPE_GEOM_ALLOC = .TRUE. IF ( I_BUG > 0 ) PRINT *, 'ALLOCATED TYPE GEOMETRY' ELSE PRINT *, 'WARNING: TYPE GEOM ALREADY ALLOCATED' N_WARN = N_WARN + 1 ! INCREMENT WARNING END IF IF ( .NOT. TYPE_SCAL_ALLOC ) THEN ALLOCATE ( H (LT_N) ) ! SCALARS ALLOCATE ( DLH (LT_PARM, LT_N) ) ALLOCATE ( D2LH (N_2_DER, LT_N) ) ALLOCATE ( H_QP (LT_N , LT_QP) ) ALLOCATE ( DLH_QP (LT_PARM, LT_N, LT_QP) ) !b IF (N_2_DER > 0) PRINT *,'NOTE: NUMBER OF 2ND DERIVATIVES =',& !b IF ( SCP_2ND_DERIV .OR. SCP_3RD_DERIV ) PRINT *, & !b 'NOTE: NUMBER OF 2ND DERIVATIVES =', N_2_DER !b IF ( SCP_3RD_DERIV ) PRINT *, & !b 'NOTE: NUMBER OF 3RD DERIVATIVES =', N_3_DER TYPE_SCAL_ALLOC = .TRUE. IF ( I_BUG > 0 ) PRINT *, 'ALLOCATED TYPE SCALARS' ELSE PRINT *, 'WARNING: TYPE SCAL ALREADY ALLOCATED' N_WARN = N_WARN + 1 ! INCREMENT WARNING END IF IF ( .NOT. TYPE_DOF_ALLOC ) THEN ALLOCATE ( V (LT_FREE) ) ! DOF ALLOCATE ( DLV (LT_PARM, LT_FREE) ) !XXX ALLOCATE ( DGV (LT_PARM, LT_FREE) ) !b LT_PARM or N_SPACE ALLOCATE ( V_QP (LT_FREE, LT_QP) ) ALLOCATE ( DLV_QP (LT_PARM, LT_FREE, LT_QP) ) TYPE_DOF_ALLOC = .TRUE. IF ( I_BUG > 0 ) PRINT *, 'ALLOCATED TYPE DOF' ELSE PRINT *, 'WARNING: TYPE NTRP ALREADY ALLOCATED' N_WARN = N_WARN + 1 ! INCREMENT WARNING END IF TYPE_NTRP_ALLOC = .TRUE. END SUBROUTINE ALLOCATE_TYPE_INTERPOLATIONS SUBROUTINE DEALLOCATE_TYPE_INTERPOLATIONS DEALLOCATE ( DLV_QP ) ! DOF DEALLOCATE ( V_QP ) !XXX DEALLOCATE ( DGV ) DEALLOCATE ( DLV ) DEALLOCATE ( V ) TYPE_DOF_ALLOC = .FALSE. IF ( I_BUG > 0 ) PRINT *, 'DEALLOCATED TYPE DOF' DEALLOCATE ( DLH_QP ) ! SCALARS DEALLOCATE ( H_QP ) DEALLOCATE ( D2LH ) DEALLOCATE ( DLH ) DEALLOCATE ( H ) TYPE_SCAL_ALLOC = .FALSE. IF ( I_BUG > 0 ) PRINT *, 'DEALLOCATED TYPE SCALARS' DEALLOCATE ( DLG_QP ) ! GEOMETRY DEALLOCATE ( G_QP ) DEALLOCATE ( DLG ) DEALLOCATE ( G ) TYPE_GEOM_ALLOC = .FALSE. IF ( I_BUG > 0 ) PRINT *, 'DEALLOCATED TYPE GEOMETRY' CALL DEALLOCATE_TYPE_QUADRATURE TYPE_NTRP_ALLOC = .FALSE. END SUBROUTINE DEALLOCATE_TYPE_INTERPOLATIONS SUBROUTINE FILL_TYPE_INTERPOLATIONS ! * * * * * * * * * * * * * * * * * * * * * * * * * * * ! FILL GEOMETRIC AND SOLUTION INTERPOLATIONS AND THEIR ! DERIVATIVES AT EACH QUADRATURE POINT ! * * * * * * * * * * * * * * * * * * * * * * * * * * * IMPLICIT NONE INTEGER :: IP IF ( LT_QP < 1 ) THEN PRINT *, 'WARNING, SKIPPING FILL_TYPE_INTERPOLATIONS' N_WARN = N_WARN + 1 ! INCREMENT WARNING ELSE ! EVALUATE GEOMETRIC INTERPOLATION FUNCTIONS DO IP = 1, LT_QP CALL GEOM_SHAPES (PT(:, IP), G_QP (:, IP)) CALL GEOM_DERIVS (PT(:, IP), DLG_QP (:, :, IP)) END DO IF ( I_BUG > 0 ) PRINT *, 'FILLED GEOMETRIC INTERPOLATIONS' ! EVALUATE SCALAR INTERPOLATION FUNCTIONS IF ( LT_GEOM == LT_N ) THEN ! SAME ARRAYS H_QP = G_QP ; DLH_QP = DLG_QP ELSE DO IP = 1, LT_QP CALL SCALAR_SHAPES (PT(:, IP), H_QP (:, IP)) CALL SCALAR_DERIVS (PT(:, IP), DLH_QP (:, :, IP)) END DO END IF ! SAME ARRAYS IF ( I_BUG > 0 ) PRINT *, 'FILLED SCALAR INTERPOLATIONS' ! EVALUATE DOF'S INTERPOLATION FUNCTIONS IF ( LT_FREE == LT_N ) THEN ! SAME ARRAYS V_QP = H_QP ; DLV_QP = DLH_QP ELSE V_QP = 0.d0 ; DLV_QP = 0.d0 print *,'NOTE: VECTOR DOF ARRAYS ARE NULL' !b IF ( N_G_DOF > 1 ) fill with copies of H, else use !b SUBROUTINE ELEM_C_N_SHAPES and SUBROUTINE ELEM_C_N_DERIVS !b or other special cases END IF ! SAME ARRAYS IF ( I_BUG > 0 ) PRINT *, 'FILLED DOF INTERPOLATIONS' END IF ! LT_QP > 0 END SUBROUTINE FILL_TYPE_INTERPOLATIONS SUBROUTINE G_AT_QP (IP) ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * ! GET LOCAL GEOMETRY INTERPOLATION AT QUADRATURE POINT ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * INTEGER :: IP ! QUADRATURE POINT NUMBER G (:) = G_QP (:, IP) END SUBROUTINE G_AT_QP FUNCTION GET_G_AT_QP (IP) RESULT (G_AT_QP) ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * ! GET LOCAL GEOMETRY INTERPOLATION AT QUADRATURE POINT ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * REAL(DP) :: G_AT_QP (LT_GEOM) INTEGER :: IP ! QUADRATURE POINT NUMBER G_AT_QP (:) = G_QP (:, IP) END FUNCTION GET_G_AT_QP FUNCTION GET_DLG_AT_QP (IP) RESULT (DLG_AT_QP) ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * ! GET LOCAL GEOMETRY DERIVATIVES AT QUADRATURE POINT ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * REAL(DP) :: DLG_AT_QP (LT_PARM, LT_GEOM) INTEGER :: IP ! QUADRATURE POINT NUMBER DLG_AT_QP (:, :) = DLG_QP (:, :, IP) END FUNCTION GET_DLG_AT_QP FUNCTION GET_H_AT_QP (IP) RESULT (H_AT_QP) ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * ! GET LOCAL SOLUTION INTERPOLATION AT QUADRATURE POINT ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * REAL(DP) :: H_AT_QP (LT_N) INTEGER :: IP ! QUADRATURE POINT NUMBER H_AT_QP (:) = H_QP (:, IP) END FUNCTION GET_H_AT_QP FUNCTION GET_DLH_AT_QP (IP) RESULT (DLH_AT_QP) ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * ! GET LOCAL SOLUTION DERIVATIVES AT QUADRATURE POINT ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * REAL(DP) :: DLH_AT_QP (LT_PARM, LT_N) INTEGER :: IP ! QUADRATURE POINT NUMBER DLH_AT_QP (:, :) = DLH_QP (:, :, IP) END FUNCTION GET_DLH_AT_QP !b XXX need GET_D2LH_AT_QP (IQ) from SCALAR_2ND_DERIVS SUBROUTINE LOCAL_COORD_AT_GP (IP, PT_IP) ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * ! GET LOCAL COORDINATES AT QUADRATURE POINT ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT (IN) :: IP REAL(DP), INTENT (OUT) :: PT_IP (LT_PARM) PT_IP (:) = PT (:, IP) END SUBROUTINE LOCAL_COORD_AT_GP SUBROUTINE WEIGHT_AT_GP (IP, WT_IP) ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * ! GET WEIGHT AT QUADRATURE POINT ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT (IN) :: IP REAL(DP), INTENT (OUT) :: WT_IP WT_IP = WT (IP) END SUBROUTINE WEIGHT_AT_GP FUNCTION GET_LOCAL_COORD_AT_GP (IP) RESULT (PT_IP) ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * ! GET LOCAL COORDINATES AT QUADRATURE POINT ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT (IN) :: IP REAL(DP) :: PT_IP (LT_PARM) PT_IP (:) = PT (:, IP) END FUNCTION GET_LOCAL_COORD_AT_GP FUNCTION GET_WEIGHT_AT_GP (IP) RESULT (WT_IP) ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * ! GET WEIGHT AT QUADRATURE POINT ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT (IN) :: IP REAL(DP) :: WT_IP WT_IP = WT (IP) END FUNCTION GET_WEIGHT_AT_GP !b SUBROUTINE ALLOCATE_ELEM_TYPE_DATA !b ALLOCATE ( L_TYPE (L_S_TOT) ) !b ALLOCATE ( LT_DATA (MAX_DAT, N_L_TYPE) ) !b TYPE_DATA_ALLOC = .TRUE. !b END SUBROUTINE ALLOCATE_ELEM_TYPE_DATA End Module Elem_Type_Data Module Sys_Properties_Data Use System_Constants IMPLICIT NONE ! flt- floating point, fix-fixed point ! b_seg-boundary segment, el-element ! misc-miscellaneous, np-nodal point REAL(DP), ALLOCATABLE :: FLT_SP (:,:), FLT_EL (:,:) REAL(DP), ALLOCATABLE :: FLT_NP (:,:), FLT_MISC (:) REAL(DP), ALLOCATABLE :: FLT_MAT (:,:), FLT_MX (:,:) REAL(DP), ALLOCATABLE :: S_FLUX (:,:,:) REAL(DP), ALLOCATABLE :: FLT_EXACT (:) INTEGER, ALLOCATABLE :: SP_FIX (:,:), LP_FIX (:,:) INTEGER, ALLOCATABLE :: NP_FIX (:,:), MISC_FIX (:) INTEGER, ALLOCATABLE :: MX_FIX (:,:) INTEGER, ALLOCATABLE :: EXACT_FIX (:) LOGICAL, ALLOCATABLE :: DELETED (:) !b LOGICAL, ALLOCATABLE :: ON_BOUNDARY (:) LOGICAL :: PROP_ELEM_ALLOC = .FALSE. LOGICAL :: PROP_EXACT_ALLOC = .FALSE. LOGICAL :: PROP_MAT_ALLOC = .FALSE. LOGICAL :: PROP_MISC_ALLOC = .FALSE. LOGICAL :: PROP_MIXED_ALLOC = .FALSE. LOGICAL :: PROP_NODE_ALLOC = .FALSE. LOGICAL :: PROP_SEGS_ALLOC = .FALSE. contains SUBROUTINE ALLOCATE_PROPERTIES ! * * * * * * * * * * * * * * * * * * * * * * * ! ALLOCATE PROPERTY ARRAYS ! * * * * * * * * * * * * * * * * * * * * * * * IF ( L_HOMO == 1 ) THEN ! homogeneous elements ALLOCATE ( LP_FIX (1, N_LP_FIX) ) ALLOCATE ( FLT_EL (1, N_LP_FLO) ) ALLOCATE ( MX_FIX (1, N_MX_FIX) ) ALLOCATE ( FLT_MX (1, N_MX_FLO) ) ELSE ALLOCATE ( LP_FIX (L_S_TOT, N_LP_FIX) ) ALLOCATE ( FLT_EL (L_S_TOT, N_LP_FLO) ) ALLOCATE ( MX_FIX (N_MIXED, N_MX_FIX) ) ALLOCATE ( FLT_MX (N_MIXED, N_MX_FLO) ) END IF ! homogeneous ALLOCATE ( DELETED (L_S_TOT)) ; DELETED = .FALSE. !b ALLOCATE ( ON_BOUNDARY (N_ELEMS)) ; ON_BOUNDARY = .FALSE. PROP_ELEM_ALLOC = .TRUE. IF ( N_HOMO == 1 ) THEN ! homogeneous nodes ALLOCATE ( NP_FIX (1, N_NP_FIX) ) ALLOCATE ( FLT_NP (1, N_NP_FLO) ) ELSE ALLOCATE ( NP_FIX (MAX_NP, N_NP_FIX)) ALLOCATE ( FLT_NP (MAX_NP, N_NP_FLO)) END IF ! homogeneous PROP_NODE_ALLOC = .TRUE. PROP_MIXED_ALLOC = .TRUE. !b ?? homogeneous el --> homogeneous seg ?? ALLOCATE ( S_FLUX (L_B_N, N_G_FLUX, N_F_SEG), & FLT_SP (N_F_SEG, N_BS_FLO), & SP_FIX (N_F_SEG, N_BS_FIX) ) PROP_SEGS_ALLOC = .TRUE. ALLOCATE ( FLT_EXACT (EXACT_FL), EXACT_FIX (EXACT_FX) ) PROP_EXACT_ALLOC = .TRUE. ALLOCATE ( FLT_MISC (MISC_FL), MISC_FIX (MISC_FX) ) PROP_MISC_ALLOC = .TRUE. ALLOCATE ( FLT_MAT (N_MAT, MAT_FLO) ) PROP_MAT_ALLOC = .TRUE. IF ( I_BUG > 0 ) PRINT *, 'ALLOCATED PROPERTIES ARRAYS' END SUBROUTINE ALLOCATE_PROPERTIES SUBROUTINE DELETE_ELEMENT (IE) INTEGER, INTENT(IN) :: IE ! Element number DELETED (IE) = .TRUE. END SUBROUTINE DELETE_ELEMENT SUBROUTINE INPUT_DELETED_ELEMENTS ! * * * * * * * * * * * * * * * * * * * * * * * ! READ DELETED ELEMENT FLAG, TERMINATE AT LAST ELEMENT ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER :: J, IE ! LOOPS LOGICAL :: TEMP ! LOOPS IF ( .NOT. PROP_ELEM_ALLOC ) THEN PRINT *, 'INPUT_DELETED_ELEMENTS: WARNING DELETED NOT ALLOCATED' N_WARN = N_WARN + 1 ! INCREMENT WARNING DO J = 1, L_S_TOT ! READ BUT IGNORE READ (N_CRD, *) IE, TEMP IF ( IE == N_ELEMS ) EXIT ! THE INPUT LOOP END DO ! OVER ELEMENTS ELSE PRINT *, 'READING INITIAL LIST OF DELETED ELEMENTS' DELETED = .FALSE. DO J = 1, L_S_TOT READ (N_CRD, *) IE, TEMP DELETED (IE) = TEMP IF ( IE == N_ELEMS ) EXIT ! THE INPUT LOOP END DO ! OVER ELEMENTS END IF END SUBROUTINE INPUT_DELETED_ELEMENTS FUNCTION IS_ELEM_DELETED (IE) RESULT (STATUS) INTEGER, INTENT(IN) :: IE ! Element number LOGICAL :: STATUS STATUS = DELETED (IE) END FUNCTION IS_ELEM_DELETED SUBROUTINE DEALLOCATE_PROPERTIES ! * * * * * * * * * * * * * * * * * * * * * * * ! DEALLOCATE PROPERTY ARRAYS, IN REVERSE ORDER ! * * * * * * * * * * * * * * * * * * * * * * * !b add status checks for each array DEALLOCATE (DELETED) DEALLOCATE (EXACT_FIX, S_FLUX, MISC_FIX, NP_FIX, SP_FIX, LP_FIX) DEALLOCATE (FLT_EXACT, FLT_MISC, FLT_NP, FLT_SP, FLT_EL) DEALLOCATE (MX_FIX, FLT_MX) !b PROP_EXACT_ALLOC = .FALSE. PROP_MISC_ALLOC = .FALSE. PROP_ELEM_ALLOC = .FALSE. PROP_NODE_ALLOC = .FALSE. PROP_SEGS_ALLOC = .FALSE. PROP_MIXED_ALLOC = .FALSE. IF ( I_BUG > 0 ) PRINT *, 'DEALLOCATED PROPERTIES ARRAYS' END SUBROUTINE DEALLOCATE_PROPERTIES SUBROUTINE SET_FIX_LP (IE, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! ASSIGN VALUE TO ELEMENT INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: IE ! Element number INTEGER, INTENT(IN) :: PROP ! Property number INTEGER, INTENT(IN) :: VALUE ! Property value INTEGER :: L_ID ! INTEGER LP_FIX (L_S_TOT, N_LP_FIX) ! in Sys_Properties_Data ! LP_FIX = SYSTEM ARRAY OF FIXED PT ELEM PROPERTIES ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! L_ID = ELEMENT NUMBER L_ID = IE ; IF ( L_HOMO == 1 ) L_ID = 1 IF ( PROP > 0 .AND. PROP <= N_LP_FIX ) THEN ! data exists LP_FIX (L_ID, PROP) = VALUE ELSE ! invalid property PRINT *,'ERROR: SET_FIX_LP, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: SET_FIX_LP, NO INPUT FOR REQUESTED PROPERTY' END IF LP_FIX (L_ID, PROP) = VALUE END SUBROUTINE SET_FIX_LP SUBROUTINE SET_FIX_MISC (PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! ASSIGN VALUE TO MISCELLANEOUS INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: PROP ! Property number INTEGER, INTENT(IN) :: VALUE ! Property value ! INTEGER MISC_FIX (MISC_FX) ! in Sys_Properties_Data ! MISC_FIX = MISCELLANEOUS INTEGER SYSTEM PROPERTIES IF ( PROP > 0 .AND. PROP <= MISC_FX ) THEN ! data exists MISC_FIX (PROP) = VALUE ELSE ! invalid property PRINT *,'ERROR: SET_FIX_MISC, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: SET_FIX_MISC, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE SET_FIX_MISC SUBROUTINE SET_FIX_MX (MX, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! ASSIGN VALUE TO MIXED_SEGMENT INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: MX ! Mixed segment number INTEGER, INTENT(IN) :: PROP ! Property number INTEGER, INTENT(IN) :: VALUE ! Property value INTEGER :: L_ID ! INTEGER MX_FIX (N_MIXED, N_MX_FIX) ! in Sys_Properties_Data ! MX_FIX = SYSTEM ARRAY OF INTEGER MIXED SEG PROPERTIES ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! L_ID = ELEMENT NUMBER L_ID = MX ; IF ( L_HOMO == 1 ) L_ID = 1 IF ( PROP > 0 .AND. PROP <= N_MX_FIX ) THEN ! data exists MX_FIX (L_ID, PROP) = VALUE ELSE ! invalid property PRINT *,'ERROR: SET_FIX_MX, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: SET_FIX_MX, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE SET_FIX_MX SUBROUTINE SET_FIX_NP (NP, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! ASSIGN VALUE TO NODAL INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: NP ! Node number INTEGER, INTENT(IN) :: PROP ! Property number INTEGER, INTENT(IN) :: VALUE ! Property value INTEGER :: N_ID ! INTEGER NP_FIX (MAX_NP, N_NP_FIX) ! in Sys_Properties_Data ! NP_FIX = SYSTEM ARRAY OF FIXED PT NODE PROPERTIES ! N_HOMO = 1, IF PROPERTIES ARE SAME AT ALL NODES ! N_ID = NODE NUMBER N_ID = NP ; IF ( N_HOMO == 1 ) N_ID = 1 IF ( PROP > 0 .AND. PROP <= N_NP_FIX ) THEN ! data exists NP_FIX (N_ID, PROP) = VALUE ELSE ! invalid property PRINT *,'ERROR: SET_FIX_NP, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: SET_FIX_NP, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE SET_FIX_NP SUBROUTINE SET_FIX_SP (NS, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! ASSIGN VALUE TO SEGMENT INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: NS ! Segment number INTEGER, INTENT(IN) :: PROP ! Property number INTEGER, INTENT(IN) :: VALUE ! Property value INTEGER :: N_ID ! INTEGER SP_FIX (N_F_SEG, N_BS_FIX) in Sys_Properties_Data ! SP_FIX = INTEGER PROPERTIES AT ALL SEGMENTS ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! N_ID = SEGMENT NUMBER N_ID = NS ; IF ( L_HOMO == 1 ) N_ID = 1 IF ( PROP > 0 .AND. PROP <= N_BS_FIX ) THEN ! data exists SP_FIX (N_ID, PROP) = VALUE ELSE ! invalid property PRINT *,'ERROR: SET_FIX_SP, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: SET_FIX_SP, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE SET_FIX_SP ! ----- rename to SET_INTEGER_... and SET_REAL_... SUBROUTINE SET_INTEGER_LP (PROP, VALUE) !b !b SUBROUTINE SET_INTEGER_LP (IE, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! ASSIGN VALUE TO ELEMENT INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * !b INTEGER, INTENT(IN) :: IE ! Element number INTEGER, INTENT(IN) :: PROP ! Property number INTEGER, INTENT(IN) :: VALUE ! Property value INTEGER :: L_ID ! INTEGER LP_FIX (L_S_TOT, N_LP_FIX) ! in Sys_Properties_Data ! LP_FIX = SYSTEM ARRAY OF FIXED PT ELEM PROPERTIES ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! L_ID = ELEMENT NUMBER !b L_ID = IE ; IF ( L_HOMO == 1 ) L_ID = 1 L_ID = GET_THIS_ELEMENT_NUMBER () ; IF ( L_HOMO == 1 ) L_ID = 1 !b IF ( PROP > 0 .AND. PROP <= N_LP_FIX ) THEN ! data exists LP_FIX (L_ID, PROP) = VALUE ELSE ! invalid property PRINT *,'ERROR: SET_INTEGER_LP, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: SET_INTEGER_LP, NO INPUT FOR REQUESTED PROPERTY' END IF LP_FIX (L_ID, PROP) = VALUE END SUBROUTINE SET_INTEGER_LP SUBROUTINE SET_INTEGER_MISC (PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! ASSIGN VALUE TO MISCELLANEOUS INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: PROP ! Property number INTEGER, INTENT(IN) :: VALUE ! Property value ! INTEGER MISC_FIX (MISC_FX) ! in Sys_Properties_Data ! MISC_FIX = MISCELLANEOUS INTEGER SYSTEM PROPERTIES IF ( PROP > 0 .AND. PROP <= MISC_FX ) THEN ! data exists MISC_FIX (PROP) = VALUE ELSE ! invalid property PRINT *,'ERROR: SET_INTEGER_MISC, NO INPUT FOR PROPERTY ', PROP STOP 'ERROR: SET_INTEGER_MISC, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE SET_INTEGER_MISC SUBROUTINE SET_INTEGER_MX (MX, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! ASSIGN VALUE TO MIXED_SEGMENT INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: MX ! Mixed segment number INTEGER, INTENT(IN) :: PROP ! Property number INTEGER, INTENT(IN) :: VALUE ! Property value INTEGER :: L_ID ! INTEGER MX_FIX (N_MIXED, N_MX_FIX) ! in Sys_Properties_Data ! MX_FIX = SYSTEM ARRAY OF INTEGER MIXED SEG PROPERTIES ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! L_ID = ELEMENT NUMBER L_ID = MX ; IF ( L_HOMO == 1 ) L_ID = 1 IF ( PROP > 0 .AND. PROP <= N_MX_FIX ) THEN ! data exists MX_FIX (L_ID, PROP) = VALUE ELSE ! invalid property PRINT *,'ERROR: SET_INTEGER_MX, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: SET_INTEGER_MX, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE SET_INTEGER_MX SUBROUTINE SET_INTEGER_NP (NP, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! ASSIGN VALUE TO NODAL INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: NP ! Node number INTEGER, INTENT(IN) :: PROP ! Property number INTEGER, INTENT(IN) :: VALUE ! Property value INTEGER :: N_ID ! INTEGER NP_FIX (MAX_NP, N_NP_FIX) ! in Sys_Properties_Data ! NP_FIX = SYSTEM ARRAY OF FIXED PT NODE PROPERTIES ! N_HOMO = 1, IF PROPERTIES ARE SAME AT ALL NODES ! N_ID = NODE NUMBER N_ID = NP ; IF ( N_HOMO == 1 ) N_ID = 1 IF ( PROP > 0 .AND. PROP <= N_NP_FIX ) THEN ! data exists NP_FIX (N_ID, PROP) = VALUE ELSE ! invalid property PRINT *,'ERROR: SET_INTEGER_NP, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: SET_INTEGER_NP, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE SET_INTEGER_NP SUBROUTINE SET_INTEGER_SP (NS, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! ASSIGN VALUE TO SEGMENT INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: NS ! Segment number INTEGER, INTENT(IN) :: PROP ! Property number INTEGER, INTENT(IN) :: VALUE ! Property value INTEGER :: N_ID ! INTEGER SP_FIX (N_F_SEG, N_BS_FIX) in Sys_Properties_Data ! SP_FIX = INTEGER PROPERTIES AT ALL SEGMENTS ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! N_ID = SEGMENT NUMBER N_ID = NS ; IF ( L_HOMO == 1 ) N_ID = 1 IF ( PROP > 0 .AND. PROP <= N_BS_FIX ) THEN ! data exists SP_FIX (N_ID, PROP) = VALUE ELSE ! invalid property PRINT *,'ERROR: SET_INTEGER_SP, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: SET_INTEGER_SP, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE SET_INTEGER_SP SUBROUTINE SET_REAL_LP (PROP, VALUE) !b !b SUBROUTINE SET_REAL_LP (IE, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! ASSIGN VALUE TO ELEMENT REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * !b INTEGER, INTENT(IN) :: IE ! Element number INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP), INTENT(IN) :: VALUE ! Property value INTEGER :: L_ID ! REAL(DP) FLT_EL (L_S_TOT, N_LP_FLO) ! in Sys_Properties_Data ! FLT_EL = SYS ARRAY OF FLOATING PT NODAL PROP ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! L_ID = ELEMENT NUMBER !b L_ID = IE ; IF ( L_HOMO == 1 ) L_ID = 1 L_ID = GET_THIS_ELEMENT_NUMBER () ; IF ( L_HOMO == 1 ) L_ID = 1 !b IF ( PROP > 0 .AND. PROP <= N_LP_FLO ) THEN ! data exists FLT_EL (L_ID, PROP) = VALUE ELSE ! invalid property PRINT *,'ERROR: SET_REAL_LP, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: SET_REAL_LP, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE SET_REAL_LP SUBROUTINE SET_REAL_MISC (PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! ASSIGN VALUE TO MISCELLANEOUS REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP), INTENT(IN) :: VALUE ! Property value ! REAL(DP) FLT_MISC (MISC_FL) ! in Sys_Properties_Data ! FLT_MISC = MISCELLANEOUS REAL PROPERTIES OF SYSTEM IF ( PROP > 0 .AND. PROP <= MISC_FL ) THEN ! data exists FLT_MISC (PROP) = VALUE ELSE ! invalid property PRINT *,'ERROR: SET_REAL_MISC, NO INPUT FOR PROPERTY ', PROP STOP 'ERROR: SET_REAL_MISC, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE SET_REAL_MISC SUBROUTINE SET_REAL_MX (MX, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! ASSIGN VALUE TO MIXED_SEGMENT REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: MX ! Mixed segment number INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP), INTENT(IN) :: VALUE ! Property value INTEGER :: L_ID ! REAL(DP) FLT_MX (N_MIXED, N_MX_FLO) ! in Sys_Properties_Data ! FLT_MX = SYS ARRAY OF REAL MIXED_SEGMENT PROPERTIES ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! L_ID = ELEMENT NUMBER L_ID = MX ; IF ( L_HOMO == 1 ) L_ID = 1 IF ( PROP > 0 .AND. PROP <= N_MX_FLO ) THEN ! data exists FLT_MX (L_ID, PROP) = VALUE ELSE ! invalid property PRINT *,'ERROR: SET_REAL_MX, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: SET_REAL_MX, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE SET_REAL_MX SUBROUTINE SET_REAL_NP (NP, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! ASSIGN VALUE TO NODAL REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: NP ! Node number INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP), INTENT(IN) :: VALUE ! Property value INTEGER :: N_ID ! REAL(DP) FLT_NP (MAX_NP, N_NP_FLO) ! in Sys_Properties_Data ! FLT_NP = REAL PROPERTIES OF SYSTEM NODES ! N_HOMO = 1, IF PROPERTIES ARE SAME AT ALL NODES ! N_ID = NODE NUMBER N_ID = NP ; IF ( N_HOMO == 1 ) N_ID = 1 IF ( PROP > 0 .AND. PROP <= N_NP_FLO ) THEN ! data exists FLT_NP (N_ID, PROP) = VALUE ELSE ! invalid property PRINT *,'ERROR: SET_REAL_NP, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: SET_REAL_NP, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE SET_REAL_NP SUBROUTINE SET_REAL_SP (NS, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! ASSIGN VALUE TO SEGMENT REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: NS ! Segment number INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP), INTENT(IN) :: VALUE ! Property value INTEGER :: N_ID ! REAL(DP) FLT_SP (N_F_SEG, N_BS_FLO) ! in Sys_Properties_Data ! FLT_SP = REAL PROPERTIES OF SYSTEM SEGMENTS ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! N_ID = SEGMENT NUMBER N_ID = NS ; IF ( L_HOMO == 1 ) N_ID = 1 IF ( PROP > 0 .AND. PROP <= N_BS_FLO ) THEN ! data exists FLT_SP (N_ID, PROP) = VALUE ELSE ! invalid property PRINT *,'ERROR: SET_REAL_SP, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: SET_REAL_SP, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE SET_REAL_SP ! renaming complete SUBROUTINE SET_FLT_LP (IE, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! ASSIGN VALUE TO ELEMENT REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: IE ! Element number INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP), INTENT(IN) :: VALUE ! Property value INTEGER :: L_ID ! REAL(DP) FLT_EL (L_S_TOT, N_LP_FLO) ! in Sys_Properties_Data ! FLT_EL = SYS ARRAY OF FLOATING PT NODAL PROP ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! L_ID = ELEMENT NUMBER L_ID = IE ; IF ( L_HOMO == 1 ) L_ID = 1 IF ( PROP > 0 .AND. PROP <= N_LP_FLO ) THEN ! data exists FLT_EL (L_ID, PROP) = VALUE ELSE ! invalid property PRINT *,'ERROR: SET_FLT_LP, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: SET_FLT_LP, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE SET_FLT_LP SUBROUTINE SET_FLT_MISC (PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! ASSIGN VALUE TO MISCELLANEOUS REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP), INTENT(IN) :: VALUE ! Property value ! REAL(DP) FLT_MISC (MISC_FL) ! in Sys_Properties_Data ! FLT_MISC = MISCELLANEOUS REAL PROPERTIES OF SYSTEM IF ( PROP > 0 .AND. PROP <= MISC_FL ) THEN ! data exists FLT_MISC (PROP) = VALUE ELSE ! invalid property PRINT *,'ERROR: SET_FLT_MISC, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: SET_FLT_MISC, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE SET_FLT_MISC SUBROUTINE SET_FLT_MX (MX, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! ASSIGN VALUE TO MIXED_SEGMENT REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: MX ! Mixed segment number INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP), INTENT(IN) :: VALUE ! Property value INTEGER :: L_ID ! REAL(DP) FLT_MX (N_MIXED, N_MX_FLO) ! in Sys_Properties_Data ! FLT_MX = SYS ARRAY OF REAL MIXED_SEGMENT PROPERTIES ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! L_ID = ELEMENT NUMBER L_ID = MX ; IF ( L_HOMO == 1 ) L_ID = 1 IF ( PROP > 0 .AND. PROP <= N_MX_FLO ) THEN ! data exists FLT_MX (L_ID, PROP) = VALUE ELSE ! invalid property PRINT *,'ERROR: SET_FLT_MX, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: SET_FLT_MX, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE SET_FLT_MX SUBROUTINE SET_FLT_NP (NP, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! ASSIGN VALUE TO NODAL REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: NP ! Node number INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP), INTENT(IN) :: VALUE ! Property value INTEGER :: N_ID ! REAL(DP) FLT_NP (MAX_NP, N_NP_FLO) ! in Sys_Properties_Data ! FLT_NP = REAL PROPERTIES OF SYSTEM NODES ! N_HOMO = 1, IF PROPERTIES ARE SAME AT ALL NODES ! N_ID = NODE NUMBER N_ID = NP ; IF ( N_HOMO == 1 ) N_ID = 1 IF ( PROP > 0 .AND. PROP <= N_NP_FLO ) THEN ! data exists FLT_NP (N_ID, PROP) = VALUE ELSE ! invalid property PRINT *,'ERROR: SET_FLT_NP, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: SET_FLT_NP, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE SET_FLT_NP SUBROUTINE SET_FLT_SP (NS, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! ASSIGN VALUE TO SEGMENT REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: NS ! Segment number INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP), INTENT(IN) :: VALUE ! Property value INTEGER :: N_ID ! REAL(DP) FLT_SP (N_F_SEG, N_BS_FLO) ! in Sys_Properties_Data ! FLT_SP = REAL PROPERTIES OF SYSTEM SEGMENTS ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! N_ID = SEGMENT NUMBER N_ID = NS ; IF ( L_HOMO == 1 ) N_ID = 1 IF ( PROP > 0 .AND. PROP <= N_BS_FLO ) THEN ! data exists FLT_SP (N_ID, PROP) = VALUE ELSE ! invalid property PRINT *,'ERROR: SET_FLT_SP, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: SET_FLT_SP, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE SET_FLT_SP FUNCTION GET_FIX_LP (IE, PROP) RESULT (VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF ELEMENT INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: IE ! Element number INTEGER, INTENT(IN) :: PROP ! Property number INTEGER :: VALUE ! Property value INTEGER :: L_ID ! Element id ! INTEGER LP_FIX (L_S_TOT, N_LP_FIX) ! in Sys_Properties_Data ! LP_FIX = SYSTEM ARRAY OF FIXED PT ELEM PROPERTIES ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! L_ID = ELEMENT NUMBER L_ID = IE ; IF ( L_HOMO == 1 ) L_ID = 1 IF ( PROP > 0 .AND. PROP <= N_LP_FIX ) THEN ! data exists VALUE = LP_FIX (L_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_FIX_LP, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: GET_FIX_LP, NO INPUT FOR REQUESTED PROPERTY' END IF END FUNCTION GET_FIX_LP FUNCTION GET_FIX_MISC (PROP) RESULT (VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF MISCELLANEOUS INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: PROP ! Property number INTEGER :: VALUE ! Property value ! INTEGER MISC_FIX (MISC_FX) ! in Sys_Properties_Data ! MISC_FIX = MISCELLANEOUS INTEGER SYSTEM PROPERTIES IF ( PROP > 0 .AND. PROP <= MISC_FX ) THEN ! data exists VALUE = MISC_FIX (PROP) ELSE ! invalid property PRINT *,'ERROR: GET_FIX_MISC, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: GET_FIX_MISC, NO INPUT FOR REQUESTED PROPERTY' END IF END FUNCTION GET_FIX_MISC FUNCTION GET_FIX_MX (MX, PROP) RESULT (VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF MIXED_SEGMENT INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: MX ! Mixed segment number INTEGER, INTENT(IN) :: PROP ! Property number INTEGER :: VALUE ! Property value INTEGER :: L_ID ! Element id ! INTEGER MX_FIX (N_MIXED, N_MX_FIX) ! in Sys_Properties_Data ! MX_FIX = SYSTEM ARRAY OF INTEGER MIXED SEG PROPERTIES ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! L_ID = ELEMENT NUMBER L_ID = MX ; IF ( L_HOMO == 1 ) L_ID = 1 IF ( PROP > 0 .AND. PROP <= N_MX_FIX ) THEN ! data exists VALUE = MX_FIX (L_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_FIX_MX, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: GET_FIX_MX, NO INPUT FOR REQUESTED PROPERTY' END IF END FUNCTION GET_FIX_MX FUNCTION GET_FIX_NP (NP, PROP) RESULT (VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF NODAL INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: NP ! Node number INTEGER, INTENT(IN) :: PROP ! Property number INTEGER :: VALUE ! Property value INTEGER :: N_ID ! Node id ! INTEGER NP_FIX (MAX_NP, N_NP_FIX) ! in Sys_Properties_Data ! NP_FIX = SYSTEM ARRAY OF FIXED PT NODE PROPERTIES ! N_HOMO = 1, IF PROPERTIES ARE SAME AT ALL NODES ! N_ID = NODE NUMBER N_ID = NP ; IF ( N_HOMO == 1 ) N_ID = 1 IF ( PROP > 0 .AND. PROP <= N_NP_FIX ) THEN ! data exists VALUE = NP_FIX (N_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_FIX_NP, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: GET_FIX_NP, NO INPUT FOR REQUESTED PROPERTY' END IF END FUNCTION GET_FIX_NP FUNCTION GET_FIX_SP (NS, PROP) RESULT (VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF SEGMENT INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: NS ! Segment number INTEGER, INTENT(IN) :: PROP ! Property number INTEGER :: VALUE ! Property value INTEGER :: N_ID ! INTEGER SP_FIX (N_F_SEG, N_BS_FIX) in Sys_Properties_Data ! SP_FIX = INTEGER PROPERTIES AT ALL SEGMENTS ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! N_ID = SEGMENT NUMBER ! NS = SEGMENT NUMBER N_ID = NS ; IF ( L_HOMO == 1 ) N_ID = 1 IF ( PROP > 0 .AND. PROP <= N_BS_FIX ) THEN ! data exists VALUE = SP_FIX (N_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_FIX_SP, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: GET_FIX_SP, NO INPUT FOR REQUESTED PROPERTY' END IF END FUNCTION GET_FIX_SP ! ----- rename to GET_INTEGER_... and GET_REAL_... FUNCTION GET_INTEGER_LP (PROP) RESULT (VALUE) !b ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF ELEMENT INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * !b INTEGER, INTENT(IN) :: IE ! Element number INTEGER, INTENT(IN) :: PROP ! Property number INTEGER :: VALUE ! Property value INTEGER :: L_ID ! Element id ! INTEGER LP_FIX (L_S_TOT, N_LP_FIX) ! in Sys_Properties_Data ! LP_FIX = SYSTEM ARRAY OF FIXED PT ELEM PROPERTIES ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! L_ID = ELEMENT NUMBER !b L_ID = IE ; IF ( L_HOMO == 1 ) L_ID = 1 L_ID = GET_THIS_ELEMENT_NUMBER () ; IF ( L_HOMO == 1 ) L_ID = 1 IF ( PROP > 0 .AND. PROP <= N_LP_FIX ) THEN ! data exists VALUE = LP_FIX (L_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_INTEGER_LP, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: GET_INTEGER_LP, NO INPUT FOR REQUESTED PROPERTY' END IF END FUNCTION GET_INTEGER_LP FUNCTION GET_INTEGER_EXACT (PROP) RESULT (VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF EXACTELLANEOUS INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: PROP ! Property number INTEGER :: VALUE ! Property value ! INTEGER EXACT_FIX (EXACT_FX) ! in Sys_Properties_Data ! EXACT_FIX = EXACTELLANEOUS INTEGER SYSTEM PROPERTIES IF ( PROP > 0 .AND. PROP <= EXACT_FX ) THEN ! data exists VALUE = EXACT_FIX (PROP) ELSE ! invalid property PRINT *,'ERROR: GET_INTEGER_EXACT, NO INPUT FOR PROPERTY ', PROP PRINT *,'CHECK exact_int IN keyword CONTROL' STOP 'ERROR: GET_INTEGER_EXACT, NO INPUT FOR REQUESTED PROPERTY' END IF END FUNCTION GET_INTEGER_EXACT FUNCTION GET_INTEGER_MISC (PROP) RESULT (VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF MISCELLANEOUS INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: PROP ! Property number INTEGER :: VALUE ! Property value ! INTEGER MISC_FIX (MISC_FX) ! in Sys_Properties_Data ! MISC_FIX = MISCELLANEOUS INTEGER SYSTEM PROPERTIES IF ( PROP > 0 .AND. PROP <= MISC_FX ) THEN ! data exists VALUE = MISC_FIX (PROP) ELSE ! invalid property PRINT *,'ERROR: GET_INTEGER_MISC, NO INPUT FOR PROPERTY ', PROP STOP 'ERROR: GET_INTEGER_MISC, NO INPUT FOR REQUESTED PROPERTY' END IF END FUNCTION GET_INTEGER_MISC FUNCTION GET_INTEGER_MX (PROP) RESULT (VALUE) !b !b FUNCTION GET_INTEGER_MX (MX, PROP) RESULT (VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF MIXED_SEGMENT INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * !b INTEGER, INTENT(IN) :: MX ! Mixed segment number INTEGER, INTENT(IN) :: PROP ! Property number INTEGER :: VALUE ! Property value INTEGER :: L_ID ! Element id ! INTEGER MX_FIX (N_MIXED, N_MX_FIX) ! in Sys_Properties_Data ! MX_FIX = SYSTEM ARRAY OF INTEGER MIXED SEG PROPERTIES ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! L_ID = ELEMENT NUMBER L_ID = IS_MIX ; IF ( L_HOMO == 1 ) L_ID = 1 !b L_ID = MX ; IF ( L_HOMO == 1 ) L_ID = 1 if (L_ID < 1) print *,'ERROR BAD MIXED SEGMENT ', L_ID !b IF ( PROP > 0 .AND. PROP <= N_MX_FIX ) THEN ! data exists VALUE = MX_FIX (L_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_INTEGER_MX, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: GET_INTEGER_MX, NO INPUT FOR REQUESTED PROPERTY' END IF END FUNCTION GET_INTEGER_MX FUNCTION GET_INTEGER_NP (NP, PROP) RESULT (VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF NODAL INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: NP ! Node number INTEGER, INTENT(IN) :: PROP ! Property number INTEGER :: VALUE ! Property value INTEGER :: N_ID ! Node id ! INTEGER NP_FIX (MAX_NP, N_NP_FIX) ! in Sys_Properties_Data ! NP_FIX = SYSTEM ARRAY OF FIXED PT NODE PROPERTIES ! N_HOMO = 1, IF PROPERTIES ARE SAME AT ALL NODES ! N_ID = NODE NUMBER N_ID = NP ; IF ( N_HOMO == 1 ) N_ID = 1 IF ( PROP > 0 .AND. PROP <= N_NP_FIX ) THEN ! data exists VALUE = NP_FIX (N_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_INTEGER_NP, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: GET_INTEGER_NP, NO INPUT FOR REQUESTED PROPERTY' END IF END FUNCTION GET_INTEGER_NP FUNCTION GET_INTEGER_SP (PROP) RESULT (VALUE) !b FUNCTION GET_INTEGER_SP (NS, PROP) RESULT (VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF SEGMENT INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * !b INTEGER, INTENT(IN) :: NS ! Segment number INTEGER, INTENT(IN) :: PROP ! Property number INTEGER :: VALUE ! Property value INTEGER :: N_ID ! INTEGER SP_FIX (N_F_SEG, N_BS_FIX) in Sys_Properties_Data ! SP_FIX = INTEGER PROPERTIES AT ALL SEGMENTS ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! N_ID = SEGMENT NUMBER ! NS = SEGMENT NUMBER N_ID = IS_SEG ; IF ( L_HOMO == 1 ) N_ID = 1 !b N_ID = NS ; IF ( L_HOMO == 1 ) N_ID = 1 IF ( PROP > 0 .AND. PROP <= N_BS_FIX ) THEN ! data exists VALUE = SP_FIX (N_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_INTEGER_SP, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: GET_INTEGER_SP, NO INPUT FOR REQUESTED PROPERTY' END IF END FUNCTION GET_INTEGER_SP FUNCTION COUNT_REAL_LP () RESULT (J) ! * * * * * * * * * * * * * * * * * * * * * * * ! COUNT NUMBER OF REAL ELEMENT PROPERTIES ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER :: J ; J = N_LP_FLO ! = el_real END FUNCTION COUNT_REAL_LP FUNCTION COUNT_INTEGER_LP () RESULT (J) ! * * * * * * * * * * * * * * * * * * * * * * * ! COUNT NUMBER OF INTEGER ELEMENT PROPERTIES ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER :: J ; J = N_LP_FIX ! = el_int END FUNCTION COUNT_INTEGER_LP FUNCTION COUNT_REAL_MAT () RESULT (J) ! * * * * * * * * * * * * * * * * * * * * * * * ! COUNT NUMBER OF REAL MATERIAL PROPERTIES ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER :: J ; J = MAT_FLO ! = mat_real END FUNCTION COUNT_REAL_MAT FUNCTION COUNT_INTEGER_MAT () RESULT (J) ! * * * * * * * * * * * * * * * * * * * * * * * ! COUNT NUMBER OF INTEGER MATERIAL PROPERTIES ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER :: J ; J = 0 ! MAT_FIX not present PRINT *,'WARNING: No integer material properties active.' N_WARN = N_WARN + 1 END FUNCTION COUNT_INTEGER_MAT FUNCTION COUNT_REAL_NP () RESULT (J) ! * * * * * * * * * * * * * * * * * * * * * * * ! COUNT NUMBER OF REAL NODAL PROPERTIES ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER :: J ; J = N_NP_FLO ! = pt_real END FUNCTION COUNT_REAL_NP FUNCTION COUNT_INTEGER_NP () RESULT (J) ! * * * * * * * * * * * * * * * * * * * * * * * ! COUNT NUMBER OF INTEGER NODAL PROPERTIES ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER :: J ; J = N_NP_FIX ! = pt_int END FUNCTION COUNT_INTEGER_NP FUNCTION COUNT_REAL_SP () RESULT (J) ! * * * * * * * * * * * * * * * * * * * * * * * ! COUNT NUMBER OF REAL SEGMENT PROPERTIES ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER :: J ; J = N_BS_FLO ! = seg_real END FUNCTION COUNT_REAL_SP FUNCTION COUNT_INTEGER_SP () RESULT (J) ! * * * * * * * * * * * * * * * * * * * * * * * ! COUNT NUMBER OF INTEGER SEGMENT PROPERTIES ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER :: J ; J = N_BS_FIX ! = seg_int END FUNCTION COUNT_INTEGER_SP FUNCTION COUNT_REAL_MX () RESULT (J) ! * * * * * * * * * * * * * * * * * * * * * * * ! COUNT NUMBER OF REAL MIXED SEG PROPERTIES ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER :: J ; J = N_MX_FLO ! = mixed_real END FUNCTION COUNT_REAL_MX FUNCTION COUNT_INTEGER_MX () RESULT (J) ! * * * * * * * * * * * * * * * * * * * * * * * ! COUNT NUMBER OF INTEGER MIXED SEG PROPERTIES ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER :: J ; J = N_MX_FIX ! = mixed_int END FUNCTION COUNT_INTEGER_MX FUNCTION GET_REAL_LP (PROP) RESULT (VALUE) !b !b FUNCTION GET_REAL_LP (IE, PROP) RESULT (VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF ELEMENT REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * !b INTEGER, INTENT(IN) :: IE ! Element number INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP) :: VALUE ! Property value INTEGER :: L_ID ! REAL(DP) FLT_EL (L_S_TOT, N_LP_FLO) ! in Sys_Properties_Data ! FLT_EL = SYS ARRAY OF FLOATING PT NODAL PROP ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! L_ID = ELEMENT NUMBER !b L_ID = IE ; IF ( L_HOMO == 1 ) L_ID = 1 L_ID = GET_THIS_ELEMENT_NUMBER () ; IF ( L_HOMO == 1 ) L_ID = 1 IF ( PROP > 0 .AND. PROP <= N_LP_FLO ) THEN ! data exists VALUE = FLT_EL (L_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_REAL_LP, NO INPUT FOR PROPERTY NUMBER ', PROP PRINT *,'CHECK el_reals IN keyword CONTROLS' STOP 'ERROR: GET_REAL_LP, NO INPUT FOR REQUESTED PROPERTY' END IF END FUNCTION GET_REAL_LP FUNCTION GET_REAL_EXACT (PROP) RESULT (VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF EXACT SOLUTION REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP) :: VALUE ! Property value ! REAL(DP) FLT_EXACT (EXACT_FL) ! in Sys_Properties_Data ! FLT_EXACT = EXACT SOLUTION REAL PROPERTIES OF SYSTEM IF ( PROP > 0 .AND. PROP <= EXACT_FL ) THEN ! data exists VALUE = FLT_EXACT (PROP) ELSE ! invalid property PRINT *,'ERROR: GET_REAL_EXACT, NO INPUT FOR PROPERTY NUMBER ', PROP PRINT *,'CHECK exact_reals IN keyword CONTROLS' STOP 'ERROR: GET_REAL_EXACT, NO INPUT FOR REQUESTED PROPERTY' END IF END FUNCTION GET_REAL_EXACT FUNCTION GET_REAL_MISC (PROP) RESULT (VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF MISCELLANEOUS REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP) :: VALUE ! Property value ! REAL(DP) FLT_MISC (MISC_FL) ! in Sys_Properties_Data ! FLT_MISC = MISCELLANEOUS REAL PROPERTIES OF SYSTEM IF ( PROP > 0 .AND. PROP <= MISC_FL ) THEN ! data exists VALUE = FLT_MISC (PROP) ELSE ! invalid property PRINT *,'ERROR: GET_REAL_MISC, NO INPUT FOR PROPERTY NUMBER ', PROP PRINT *,'CHECK reals IN keyword CONTROLS' STOP 'ERROR: GET_REAL_MISC, NO INPUT FOR REQUESTED PROPERTY' END IF END FUNCTION GET_REAL_MISC FUNCTION GET_REAL_MX (PROP) RESULT (VALUE) !b !b FUNCTION GET_REAL_MX (MX, PROP) RESULT (VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF MIXED_SEGMENT REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * !b INTEGER, INTENT(IN) :: MX ! Mixed segment number INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP) :: VALUE ! Property value INTEGER :: L_ID ! INTEGER :: FLT_MX (N_MIXED, N_MX_FLO) ! in Sys_Properties_Data ! FLT_MX = ARRAY OF REAL MIXED_SEGMENT PROPERTIES ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! L_ID = ELEMENT NUMBER L_ID = IS_MIX ; IF ( L_HOMO == 1 ) L_ID = 1 !b if (L_ID < 1) print *,'ERROR BAD MIXED SEGMENT ', L_ID !b !B L_ID = MX ; IF ( L_HOMO == 1 ) L_ID = 1 IF ( PROP > 0 .AND. PROP <= N_MX_FLO ) THEN ! data exists VALUE = FLT_MX (L_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_REAL_MX, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: GET_REAL_MX, NO INPUT FOR REQUESTED PROPERTY' END IF END FUNCTION GET_REAL_MX FUNCTION GET_REAL_NP (NP, PROP) RESULT (VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF NODAL REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: NP ! Node number INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP) :: VALUE ! Property value INTEGER :: N_ID ! REAL(DP) FLT_NP (MAX_NP, N_NP_FLO) ! in Sys_Properties_Data ! FLT_NP = REAL PROPERTIES OF SYSTEM NODES ! N_HOMO = 1, IF PROPERTIES ARE SAME AT ALL NODES ! N_ID = NODE NUMBER N_ID = NP ; IF ( N_HOMO == 1 ) N_ID = 1 IF ( PROP > 0 .AND. PROP <= N_NP_FLO ) THEN ! data exists VALUE = FLT_NP (N_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_REAL_NP, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: GET_REAL_NP, NO INPUT FOR REQUESTED PROPERTY' END IF END FUNCTION GET_REAL_NP FUNCTION GET_REAL_SP (PROP) RESULT (VALUE) !b FUNCTION GET_REAL_SP (NS, PROP) RESULT (VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF SEGMENT REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * !b INTEGER, INTENT(IN) :: NS ! Segment number INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP) :: VALUE ! Property value INTEGER :: N_ID ! REAL(DP) FLT_SP (N_F_SEG, N_BS_FLO) ! in Sys_Properties_Data ! FLT_SP = REAL PROPERTIES OF SYSTEM SEGMENTS ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! N_ID = SEGMENT NUMBER N_ID = IS_SEG ; IF ( L_HOMO == 1 ) N_ID = 1 !b N_ID = NS ; IF ( L_HOMO == 1 ) N_ID = 1 IF ( PROP > 0 .AND. PROP <= N_BS_FLO ) THEN ! data exists VALUE = FLT_SP (N_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_REAL_SP, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: GET_REAL_SP, NO INPUT FOR REQUESTED PROPERTY' END IF END FUNCTION GET_REAL_SP ! end rename FUNCTION GET_FLT_LP (IE, PROP) RESULT (VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF ELEMENT REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: IE ! Element number INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP) :: VALUE ! Property value INTEGER :: L_ID ! REAL(DP) FLT_EL (L_S_TOT, N_LP_FLO) ! in Sys_Properties_Data ! FLT_EL = SYS ARRAY OF FLOATING PT NODAL PROP ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! L_ID = ELEMENT NUMBER L_ID = IE ; IF ( L_HOMO == 1 ) L_ID = 1 IF ( PROP > 0 .AND. PROP <= N_LP_FLO ) THEN ! data exists VALUE = FLT_EL (L_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_FLT_LP, NO INPUT FOR PROPERTY NUMBER ', PROP PRINT *,'CHECK el_reals IN keyword CONTROLS' STOP 'ERROR: GET_FLT_LP, NO INPUT FOR REQUESTED PROPERTY' END IF END FUNCTION GET_FLT_LP FUNCTION GET_FLT_MISC (PROP) RESULT (VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF MISCELLANEOUS REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP) :: VALUE ! Property value ! REAL(DP) FLT_MISC (MISC_FL) ! in Sys_Properties_Data ! FLT_MISC = MISCELLANEOUS REAL PROPERTIES OF SYSTEM IF ( PROP > 0 .AND. PROP <= MISC_FL ) THEN ! data exists VALUE = FLT_MISC (PROP) ELSE ! invalid property PRINT *,'ERROR: GET_FLT_MISC, NO INPUT FOR PROPERTY NUMBER ', PROP PRINT *,'CHECK reals IN keyword CONTROLS' STOP 'ERROR: GET_FLT_MISC, NO INPUT FOR REQUESTED PROPERTY' END IF END FUNCTION GET_FLT_MISC FUNCTION GET_FLT_MX (MX, PROP) RESULT (VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF MIXED_SEGMENT REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: MX ! Mixed segment number INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP) :: VALUE ! Property value INTEGER :: L_ID ! INTEGER :: FLT_MX (N_MIXED, N_MX_FLO) ! in Sys_Properties_Data ! FLT_MX = ARRAY OF REAL MIXED_SEGMENT PROPERTIES ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! L_ID = ELEMENT NUMBER L_ID = MX ; IF ( L_HOMO == 1 ) L_ID = 1 IF ( PROP > 0 .AND. PROP <= N_MX_FLO ) THEN ! data exists VALUE = FLT_MX (L_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_FLT_MX, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: GET_FLT_MX, NO INPUT FOR REQUESTED PROPERTY' END IF END FUNCTION GET_FLT_MX FUNCTION GET_FLT_NP (NP, PROP) RESULT (VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF NODAL REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: NP ! Node number INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP) :: VALUE ! Property value INTEGER :: N_ID ! REAL(DP) FLT_NP (MAX_NP, N_NP_FLO) ! in Sys_Properties_Data ! FLT_NP = REAL PROPERTIES OF SYSTEM NODES ! N_HOMO = 1, IF PROPERTIES ARE SAME AT ALL NODES ! N_ID = NODE NUMBER N_ID = NP ; IF ( N_HOMO == 1 ) N_ID = 1 IF ( PROP > 0 .AND. PROP <= N_NP_FLO ) THEN ! data exists VALUE = FLT_NP (N_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_FLT_NP, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: GET_FLT_NP, NO INPUT FOR REQUESTED PROPERTY' END IF END FUNCTION GET_FLT_NP FUNCTION GET_FLT_SP (NS, PROP) RESULT (VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF SEGMENT REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: NS ! Segment number INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP) :: VALUE ! Property value INTEGER :: N_ID ! REAL(DP) FLT_SP (N_F_SEG, N_BS_FLO) ! in Sys_Properties_Data ! FLT_SP = REAL PROPERTIES OF SYSTEM SEGMENTS ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! N_ID = SEGMENT NUMBER N_ID = NS ; IF ( L_HOMO == 1 ) N_ID = 1 IF ( PROP > 0 .AND. PROP <= N_BS_FLO ) THEN ! data exists VALUE = FLT_SP (N_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_FLT_SP, NO INPUT FOR PROPERTY NUMBER ', PROP STOP 'ERROR: GET_FLT_SP, NO INPUT FOR REQUESTED PROPERTY' END IF END FUNCTION GET_FLT_SP SUBROUTINE GET_PROP_FIX_LP (IE, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF ELEMENT INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: IE ! Element number INTEGER, INTENT(IN) :: PROP ! Property number INTEGER, INTENT(OUT):: VALUE ! Property value INTEGER :: L_ID ! INTEGER LP_FIX (L_S_TOT, N_LP_FIX) ! in Sys_Properties_Data ! LP_FIX = SYSTEM ARRAY OF FIXED PT ELEM PROPERTIES ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! L_ID = ELEMENT NUMBER L_ID = IE ; IF ( L_HOMO == 1 ) L_ID = 1 VALUE = LP_FIX (L_ID, PROP) END SUBROUTINE GET_PROP_FIX_LP SUBROUTINE GET_PROP_FIX_MISC (PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF MISCELLANEOUS INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: PROP ! Property number INTEGER, INTENT(OUT):: VALUE ! Property value ! INTEGER MISC_FIX (MISC_FX) ! in Sys_Properties_Data ! MISC_FIX = MISCELLANEOUS INTEGER SYSTEM PROPERTIES IF ( PROP > 0 .AND. PROP <= MISC_FX ) THEN ! data exists VALUE = MISC_FIX (PROP) ELSE ! invalid property PRINT *,'ERROR: GET_PROP_FIX_MISC, NO INPUT FOR PROPERTY ', PROP STOP 'ERROR: GET_PROP_FIX_MISC, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE GET_PROP_FIX_MISC SUBROUTINE GET_PROP_FIX_MX (MX, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF MIXED_SEGMENT INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: MX ! Mixed segment number INTEGER, INTENT(IN) :: PROP ! Property number INTEGER, INTENT(OUT):: VALUE ! Property value INTEGER :: L_ID ! INTEGER :: MX_FIX (N_MIXED, N_MX_FIX) ! in Sys_Properties_Data ! MX_FIX = ARRAY OF INTEGER MIXED_SEG PROPERTIES ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! L_ID = ELEMENT NUMBER L_ID = MX ; IF ( L_HOMO == 1 ) L_ID = 1 IF ( PROP > 0 .AND. PROP <= N_MX_FIX ) THEN ! data exists VALUE = MX_FIX (L_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_PROP_FIX_MX, NO INPUT FOR PROPERTY ', PROP STOP 'ERROR: GET_PROP_FIX_MX, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE GET_PROP_FIX_MX SUBROUTINE GET_PROP_FIX_NP (NP, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF NODAL INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: NP ! Node number INTEGER, INTENT(IN) :: PROP ! Property number INTEGER, INTENT(OUT):: VALUE ! Property value INTEGER :: N_ID ! INTEGER NP_FIX (MAX_NP, N_NP_FIX) ! in Sys_Properties_Data ! NP_FIX = SYSTEM ARRAY OF FIXED PT NODE PROPERTIES ! N_HOMO = 1, IF PROPERTIES ARE SAME AT ALL NODES ! N_ID = NODE NUMBER N_ID = NP ; IF ( N_HOMO == 1 ) N_ID = 1 IF ( PROP > 0 .AND. PROP <= N_NP_FIX ) THEN ! data exists VALUE = NP_FIX (N_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_PROP_FIX_NP, NO INPUT FOR PROPERTY ', PROP STOP 'ERROR: GET_PROP_FIX_NP, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE GET_PROP_FIX_NP SUBROUTINE GET_PROP_FIX_SP (NS, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF SEGMENT INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: NS ! Segment number INTEGER, INTENT(IN) :: PROP ! Property number INTEGER, INTENT(OUT):: VALUE ! Property value INTEGER :: N_ID ! INTEGER SP_FIX (N_F_SEG, N_BS_FIX) in Sys_Properties_Data ! SP_FIX = INTEGER PROPERTIES AT ALL SEGMENTS ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! N_ID = SEGMENT NUMBER N_ID = NS ; IF ( L_HOMO == 1 ) N_ID = 1 IF ( PROP > 0 .AND. PROP <= N_BS_FIX ) THEN ! data exists VALUE = SP_FIX (N_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_PROP_FIX_SP, NO INPUT FOR PROPERTY ', PROP STOP 'ERROR: GET_PROP_FIX_SP, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE GET_PROP_FIX_SP ! ----- rename to GET_PROP_INTEGER_... and GET_PROP_REAL_... SUBROUTINE GET_PROP_INTEGER_LP (IE, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF ELEMENT INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: IE ! Element number INTEGER, INTENT(IN) :: PROP ! Property number INTEGER, INTENT(OUT):: VALUE ! Property value INTEGER :: L_ID ! INTEGER LP_FIX (L_S_TOT, N_LP_FIX) ! in Sys_Properties_Data ! LP_FIX = SYSTEM ARRAY OF FIXED PT ELEM PROPERTIES ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! L_ID = ELEMENT NUMBER L_ID = IE ; IF ( L_HOMO == 1 ) L_ID = 1 VALUE = LP_FIX (L_ID, PROP) END SUBROUTINE GET_PROP_INTEGER_LP SUBROUTINE GET_PROP_INTEGER_EXACT (PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF EXACTELLANEOUS INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: PROP ! Property number INTEGER, INTENT(OUT):: VALUE ! Property value ! INTEGER EXACT_FIX (EXACT_FX) ! in Sys_Properties_Data ! EXACT_FIX = EXACTELLANEOUS INTEGER SYSTEM PROPERTIES IF ( PROP > 0 .AND. PROP <= EXACT_FX ) THEN ! data exists VALUE = EXACT_FIX (PROP) ELSE ! invalid property PRINT *,'ERROR: GET_PROP_INTEGER_EXACT, NO INPUT FOR PROPERTY', PROP PRINT *,'CHECK exact_int IN keyword CONTROLS' STOP 'ERROR: GET_PROP_INTEGER_EXACT, NO INPUT FOR PROPERTY' END IF END SUBROUTINE GET_PROP_INTEGER_EXACT SUBROUTINE GET_PROP_INTEGER_MISC (PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF MISCELLANEOUS INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: PROP ! Property number INTEGER, INTENT(OUT):: VALUE ! Property value ! INTEGER MISC_FIX (MISC_FX) ! in Sys_Properties_Data ! MISC_FIX = MISCELLANEOUS INTEGER SYSTEM PROPERTIES IF ( PROP > 0 .AND. PROP <= MISC_FX ) THEN ! data exists VALUE = MISC_FIX (PROP) ELSE ! invalid property PRINT *,'ERROR: GET_PROP_INTEGER_MISC, NO INPUT FOR PROPERTY ', PROP STOP 'ERROR: GET_PROP_INTEGER_MISC, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE GET_PROP_INTEGER_MISC SUBROUTINE GET_PROP_INTEGER_MX (MX, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF MIXED_SEGMENT INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: MX ! Mixed segment number INTEGER, INTENT(IN) :: PROP ! Property number INTEGER, INTENT(OUT):: VALUE ! Property value INTEGER :: L_ID ! INTEGER :: MX_FIX (N_MIXED, N_MX_FIX) ! in Sys_Properties_Data ! MX_FIX = ARRAY OF INTEGER MIXED_SEG PROPERTIES ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! L_ID = ELEMENT NUMBER L_ID = MX ; IF ( L_HOMO == 1 ) L_ID = 1 IF ( PROP > 0 .AND. PROP <= N_MX_FIX ) THEN ! data exists VALUE = MX_FIX (L_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_PROP_INTEGER_MX, NO INPUT FOR PROPERTY ', PROP STOP 'ERROR: GET_PROP_INTEGER_MX, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE GET_PROP_INTEGER_MX SUBROUTINE GET_PROP_INTEGER_NP (NP, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF NODAL INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: NP ! Node number INTEGER, INTENT(IN) :: PROP ! Property number INTEGER, INTENT(OUT):: VALUE ! Property value INTEGER :: N_ID ! INTEGER NP_FIX (MAX_NP, N_NP_FIX) ! in Sys_Properties_Data ! NP_FIX = SYSTEM ARRAY OF FIXED PT NODE PROPERTIES ! N_HOMO = 1, IF PROPERTIES ARE SAME AT ALL NODES ! N_ID = NODE NUMBER N_ID = NP ; IF ( N_HOMO == 1 ) N_ID = 1 IF ( PROP > 0 .AND. PROP <= N_NP_FIX ) THEN ! data exists VALUE = NP_FIX (N_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_PROP_INTEGER_NP, NO INPUT FOR PROPERTY ', PROP STOP 'ERROR: GET_PROP_INTEGER_NP, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE GET_PROP_INTEGER_NP SUBROUTINE GET_PROP_INTEGER_SP (NS, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF SEGMENT INTEGER PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: NS ! Segment number INTEGER, INTENT(IN) :: PROP ! Property number INTEGER, INTENT(OUT):: VALUE ! Property value INTEGER :: N_ID ! INTEGER SP_FIX (N_F_SEG, N_BS_FIX) in Sys_Properties_Data ! SP_FIX = INTEGER PROPERTIES AT ALL SEGMENTS ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! N_ID = SEGMENT NUMBER N_ID = NS ; IF ( L_HOMO == 1 ) N_ID = 1 IF ( PROP > 0 .AND. PROP <= N_BS_FIX ) THEN ! data exists VALUE = SP_FIX (N_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_PROP_INTEGER_SP, NO INPUT FOR PROPERTY ', PROP STOP 'ERROR: GET_PROP_INTEGER_SP, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE GET_PROP_INTEGER_SP SUBROUTINE GET_PROP_REAL_LP (IE, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF ELEMENT REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: IE ! Element number INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP), INTENT(OUT):: VALUE ! Property value INTEGER :: L_ID ! REAL(DP) FLT_EL (L_S_TOT, N_LP_FLO) ! in Sys_Properties_Data ! FLT_EL = SYS ARRAY OF FLOATING PT NODAL PROP ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! L_ID = ELEMENT NUMBER L_ID = IE ; IF ( L_HOMO == 1 ) L_ID = 1 IF ( PROP > 0 .AND. PROP <= N_LP_FLO ) THEN ! data exists VALUE = FLT_EL (L_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_PROP_REAL_LP, NO INPUT FOR PROPERTY ', PROP STOP 'ERROR: GET_PROP_REAL_LP, NO INPUT FOR REQUESTED PROPERTY' END IF VALUE = FLT_EL (L_ID, PROP) END SUBROUTINE GET_PROP_REAL_LP SUBROUTINE GET_PROP_REAL_EXACT (PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF EXACT SOLUTION REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP), INTENT(OUT):: VALUE ! Property value ! REAL(DP) FLT_EXACT (EXACT_FL) ! in Sys_Properties_Data ! FLT_EXACT = EXACT SOLUTION REAL PROPERTIES OF SYSTEM IF ( PROP > 0 .AND. PROP <= EXACT_FL ) THEN ! data exists VALUE = FLT_EXACT (PROP) ELSE ! invalid property PRINT *,'ERROR: GET_PROP_REAL_EXACT, NO INPUT FOR PROPERTY ', PROP PRINT *,'CHECK exact_reals IN keyword CONTROLS' STOP 'ERROR: GET_PROP_REAL_EXACT, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE GET_PROP_REAL_EXACT SUBROUTINE GET_PROP_REAL_MISC (PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF MISCELLANEOUS REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP), INTENT(OUT):: VALUE ! Property value ! REAL(DP) FLT_MISC (MISC_FL) ! in Sys_Properties_Data ! FLT_MISC = MISCELLANEOUS REAL PROPERTIES OF SYSTEM IF ( PROP > 0 .AND. PROP <= MISC_FL ) THEN ! data exists VALUE = FLT_MISC (PROP) ELSE ! invalid property PRINT *,'ERROR: GET_PROP_REAL_MISC, NO INPUT FOR PROPERTY ', PROP STOP 'ERROR: GET_PROP_REAL_MISC, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE GET_PROP_REAL_MISC SUBROUTINE GET_PROP_REAL_MX (MX, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF MIXED_SEGMENT REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: MX ! Mixed segment number INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP), INTENT(OUT):: VALUE ! Property value INTEGER :: L_ID ! INTEGER :: FLT_MX (N_MIXED, N_MX_FLO) ! in Sys_Properties_Data ! FLT_MX = ARRAY OF REAL MIXED_SEGMENT PROPERTIES ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! L_ID = ELEMENT NUMBER L_ID = MX ; IF ( L_HOMO == 1 ) L_ID = 1 IF ( PROP > 0 .AND. PROP <= N_MX_FLO ) THEN ! data exists VALUE = FLT_MX (L_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_PROP_REAL_MX, NO INPUT FOR PROPERTY ', PROP STOP 'ERROR: GET_PROP_REAL_MX, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE GET_PROP_REAL_MX SUBROUTINE GET_PROP_REAL_NP (NP, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF NODAL REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: NP ! Node number INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP), INTENT(OUT):: VALUE ! Property value INTEGER :: N_ID ! REAL(DP) FLT_NP (MAX_NP, N_NP_FLO) ! in Sys_Properties_Data ! FLT_NP = REAL PROPERTIES OF SYSTEM NODES ! N_HOMO = 1, IF PROPERTIES ARE SAME AT ALL NODES ! N_ID = NODE NUMBER N_ID = NP ; IF ( N_HOMO == 1 ) N_ID = 1 IF ( PROP > 0 .AND. PROP <= N_NP_FLO ) THEN ! data exists VALUE = FLT_NP (N_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_PROP_REAL_NP, NO INPUT FOR PROPERTY ', PROP STOP 'ERROR: GET_PROP_REAL_NP, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE GET_PROP_REAL_NP SUBROUTINE GET_PROP_REAL_SP (NS, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF SEGMENT REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: NS ! Segment number INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP), INTENT(OUT):: VALUE ! Property value INTEGER :: N_ID ! REAL(DP) FLT_SP (N_F_SEG, N_BS_FLO) ! in Sys_Properties_Data ! FLT_SP = REAL PROPERTIES OF SYSTEM SEGMENTS ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! N_ID = SEGMENT NUMBER N_ID = NS ; IF ( L_HOMO == 1 ) N_ID = 1 IF ( PROP > 0 .AND. PROP <= N_BS_FLO ) THEN ! data exists VALUE = FLT_SP (N_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_PROP_REAL_SP, NO INPUT FOR PROPERTY ', PROP STOP 'ERROR: GET_PROP_REAL_SP, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE GET_PROP_REAL_SP ! end rename SUBROUTINE GET_PROP_FLT_LP (IE, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF ELEMENT REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: IE ! Element number INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP), INTENT(OUT):: VALUE ! Property value INTEGER :: L_ID ! REAL(DP) FLT_EL (L_S_TOT, N_LP_FLO) ! in Sys_Properties_Data ! FLT_EL = SYS ARRAY OF FLOATING PT NODAL PROP ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! L_ID = ELEMENT NUMBER L_ID = IE ; IF ( L_HOMO == 1 ) L_ID = 1 IF ( PROP > 0 .AND. PROP <= N_LP_FLO ) THEN ! data exists VALUE = FLT_EL (L_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_PROP_FLT_LP, NO INPUT FOR PROPERTY ', PROP STOP 'ERROR: GET_PROP_FLT_LP, NO INPUT FOR REQUESTED PROPERTY' END IF VALUE = FLT_EL (L_ID, PROP) END SUBROUTINE GET_PROP_FLT_LP SUBROUTINE GET_PROP_FLT_MISC (PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF MISCELLANEOUS REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP), INTENT(OUT):: VALUE ! Property value ! REAL(DP) FLT_MISC (MISC_FL) ! in Sys_Properties_Data ! FLT_MISC = MISCELLANEOUS REAL PROPERTIES OF SYSTEM IF ( PROP > 0 .AND. PROP <= MISC_FL ) THEN ! data exists VALUE = FLT_MISC (PROP) ELSE ! invalid property PRINT *,'ERROR: GET_PROP_FLT_MISC, NO INPUT FOR PROPERTY ', PROP STOP 'ERROR: GET_PROP_FLT_MISC, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE GET_PROP_FLT_MISC SUBROUTINE GET_PROP_FLT_MX (MX, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF MIXED_SEGMENT REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: MX ! Mixed segment number INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP), INTENT(OUT):: VALUE ! Property value INTEGER :: L_ID ! INTEGER :: FLT_MX (N_MIXED, N_MX_FLO) ! in Sys_Properties_Data ! FLT_MX = ARRAY OF REAL MIXED_SEGMENT PROPERTIES ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! L_ID = ELEMENT NUMBER L_ID = MX ; IF ( L_HOMO == 1 ) L_ID = 1 IF ( PROP > 0 .AND. PROP <= N_MX_FLO ) THEN ! data exists VALUE = FLT_MX (L_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_PROP_FLT_MX, NO INPUT FOR PROPERTY ', PROP STOP 'ERROR: GET_PROP_FLT_MX, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE GET_PROP_FLT_MX SUBROUTINE GET_PROP_FLT_NP (NP, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF NODAL REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: NP ! Node number INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP), INTENT(OUT):: VALUE ! Property value INTEGER :: N_ID ! REAL(DP) FLT_NP (MAX_NP, N_NP_FLO) ! in Sys_Properties_Data ! FLT_NP = REAL PROPERTIES OF SYSTEM NODES ! N_HOMO = 1, IF PROPERTIES ARE SAME AT ALL NODES ! N_ID = NODE NUMBER N_ID = NP ; IF ( N_HOMO == 1 ) N_ID = 1 IF ( PROP > 0 .AND. PROP <= N_NP_FLO ) THEN ! data exists VALUE = FLT_NP (N_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_PROP_FLT_NP, NO INPUT FOR PROPERTY ', PROP STOP 'ERROR: GET_PROP_FLT_NP, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE GET_PROP_FLT_NP SUBROUTINE GET_PROP_FLT_SP (NS, PROP, VALUE) ! * * * * * * * * * * * * * * * * * * * * * * * ! GET VALUE OF SEGMENT REAL PROPERTY ! * * * * * * * * * * * * * * * * * * * * * * * INTEGER, INTENT(IN) :: NS ! Segment number INTEGER, INTENT(IN) :: PROP ! Property number REAL(DP), INTENT(OUT):: VALUE ! Property value INTEGER :: N_ID ! REAL(DP) FLT_SP (N_F_SEG, N_BS_FLO) ! in Sys_Properties_Data ! FLT_SP = REAL PROPERTIES OF SYSTEM SEGMENTS ! L_HOMO = 1, IF PROPERTIES ARE SAME IN ALL ELEMENTS ! N_ID = SEGMENT NUMBER N_ID = NS ; IF ( L_HOMO == 1 ) N_ID = 1 IF ( PROP > 0 .AND. PROP <= N_BS_FLO ) THEN ! data exists VALUE = FLT_SP (N_ID, PROP) ELSE ! invalid property PRINT *,'ERROR: GET_PROP_FLT_SP, NO INPUT FOR PROPERTY ', PROP STOP 'ERROR: GET_PROP_FLT_SP, NO INPUT FOR REQUESTED PROPERTY' END IF END SUBROUTINE GET_PROP_FLT_SP End Module Sys_Properties_Data MODULE INTERFACE_HEADER Use Precision_Module ! NOTE: Functions returning an array require an interface INTERFACE FUNCTION COPY_DGH_INTO_B_MATRIX (DGH) RESULT (B) Use System_Constants Use Elem_Type_Data ! for LT_N, LT_FREE REAL (DP), INTENT(IN) :: DGH (N_SPACE, LT_N) REAL (DP) :: B (N_R_B, LT_FREE) END FUNCTION COPY_DGH_INTO_B_MATRIX FUNCTION GEOMETRIC_JACOBIAN () RESULT (AJ) Use System_Constants ! For N_SPACE Use Elem_Type_Data ! For LT_GEOM, LT_PARM REAL(DP) :: AJ (LT_PARM, LT_PARM) ! Out END FUNCTION GEOMETRIC_JACOBIAN FUNCTION GET_APPLICATION_B_MATRIX (DGH, XYZ) RESULT (B) Use System_Constants Use Elem_Type_Data ! for H (LT_N) IF NEEDED REAL (DP), INTENT(IN) :: DGH (N_SPACE, LT_N) REAL (DP), INTENT(IN) :: XYZ (N_SPACE) REAL (DP) :: B (N_R_B, LT_FREE) END FUNCTION GET_APPLICATION_B_MATRIX FUNCTION GET_APPLICATION_E_MATRIX (IE, XYZ) RESULT (E) Use System_Constants Use Sys_Properties_Data INTEGER, INTENT(IN) :: IE REAL (DP), INTENT(IN) :: XYZ (N_SPACE) REAL (DP) :: E (N_R_B, N_R_B) END FUNCTION GET_APPLICATION_E_MATRIX FUNCTION GET_PT_DOF_INDEX (I_PT, J_PARM) RESULT (INDEX) Use System_Constants ! for N_G_DOF INTEGER, INTENT(IN) :: I_PT, J_PARM INTEGER :: INDEX END FUNCTION GET_PT_DOF_INDEX FUNCTION GET_ELEM_DOF (DD) RESULT (D_LT) Use System_Constants ! for N_D_FREE Use Elem_Type_Data ! for INDEX (LT_FREE) REAL(DP), INTENT(IN) :: DD (N_D_FRE) REAL(DP) :: D_LT (LT_FREE) END FUNCTION GET_ELEM_DOF FUNCTION GET_ELEM_INDEX (LT_N, ELEM_NODES) RESULT (INDEX) Use System_Constants INTEGER, INTENT(IN) :: LT_N, ELEM_NODES (LT_N) INTEGER :: INDEX (LT_N*N_G_DOF) END FUNCTION GET_ELEM_INDEX FUNCTION GET_ELEM_NODES (L_ID, LT_N, NODES) RESULT (ELEM_NODES) Use System_Constants INTEGER, INTENT(IN) :: NODES (L_S_TOT, NOD_PER_EL), L_ID, LT_N INTEGER :: ELEM_NODES (LT_N) END FUNCTION GET_ELEM_NODES FUNCTION GET_INDEX_AT_PT (I_PT) RESULT (INDEX) Use System_Constants ! for N_G_DOF INTEGER, INTENT(IN) :: I_PT INTEGER :: INDEX (N_G_DOF) END FUNCTION GET_INDEX_AT_PT SUBROUTINE GET_NUMB_QUADRATURES (PT_N, WT_N, N_GP) Use Precision_Module Use Elem_Type_Data INTEGER, INTENT(IN) :: N_GP REAL(DP), INTENT(OUT) :: PT_N (LT_PARM, N_GP), WT_N (N_GP) END SUBROUTINE GET_NUMB_QUADRATURES FUNCTION GET_PT_COORD (I_PT, X) RESULT (COORD_PT) Use System_Constants INTEGER, INTENT(IN) :: I_PT REAL(DP), INTENT(IN) :: X (MAX_NP, N_SPACE) REAL(DP) :: COORD_PT (N_SPACE) END FUNCTION GET_PT_COORD FUNCTION GET_REAL_IDENTITY (N) RESULT (EYE) Use System_Constants INTEGER, INTENT(IN) :: N ! SIZE OF MATRIX REAL(DP) :: EYE (N, N) ! IDENTITY MATRIX END FUNCTION GET_REAL_IDENTITY FUNCTION GET_SKY_DIAG (I_DOF_HI) RESULT (I_DIAG) Use System_Constants INTEGER, INTENT(IN) :: I_DOF_HI (N_D_FRE) INTEGER :: I_DIAG (N_D_FRE) END FUNCTION GET_SKY_DIAG FUNCTION LOCATION_IN_SKY (I, J, I_DIAG) RESULT (IN_SKY) Use system_constants ! N_D_FRE INTEGER, INTENT(IN) :: I, J ! full subscripts INTEGER, INTENT(IN) :: I_DIAG (N_D_FRE) ! diagonal locations INTEGER :: IN_SKY ! location in skyline END FUNCTION LOCATION_IN_SKY FUNCTION LAST_COL_AT_SKY_ROW (ROW, I_DIAG) RESULT (COLUMN) Use system_constants ! N_D_FRE INTEGER, INTENT(IN) :: ROW ! row (or column) number INTEGER, INTENT(IN) :: I_DIAG (N_D_FRE) ! skyline diagonals INTEGER :: COLUMN ! first column ( or row) END FUNCTION LAST_COL_AT_SKY_ROW FUNCTION FIRST_COL_AT_SKY_ROW (ROW, I_DIAG) RESULT (COLUMN) Use system_constants ! N_D_FRE INTEGER, INTENT(IN) :: ROW ! row (or column) number INTEGER, INTENT(IN) :: I_DIAG (N_D_FRE) ! skyline diagonals INTEGER :: COLUMN ! first column ( or row) END FUNCTION FIRST_COL_AT_SKY_ROW FUNCTION OUTER_PRODUCT (A, B) RESULT (C) Use Precision_Module REAL(DP), INTENT(IN) :: A (:), B (:) REAL(DP) :: C (SIZE(A), SIZE(B)) END FUNCTION OUTER_PRODUCT SUBROUTINE SYMMETRIC_INVERSE (A, A_INV, I_ERROR) Use Precision_Module REAL (DP), INTENT (IN) :: A (:, :) REAL (DP), INTENT (OUT) :: A_INV (SIZE(A, 1), SIZE(A, 1)) INTEGER, INTENT (OUT) :: I_ERROR END SUBROUTINE SYMMETRIC_INVERSE SUBROUTINE INVERT_SMALL_MATRIX (A, A_INV, I_ERROR) Use Precision_Module REAL (DP), INTENT (IN) :: A (:, :) REAL (DP), INTENT (OUT) :: A_INV (SIZE(A, 1), SIZE(A, 1)) INTEGER, INTENT (OUT) :: I_ERROR END SUBROUTINE INVERT_SMALL_MATRIX END INTERFACE END MODULE INTERFACE_HEADER MODULE GEOMETRIC_PROPERTIES Use System_Constants Use Elem_Type_Data Use Interface_Header IMPLICIT NONE REAL(DP) :: VOLUME REAL(DP), ALLOCATABLE :: FIRST_MOMENTS (:) REAL(DP), ALLOCATABLE :: CENTROID (:) REAL(DP), ALLOCATABLE :: SECOND_MOMENTS (:, :) REAL(DP), ALLOCATABLE :: CG_SECOND_MOM (:, :) REAL(DP) :: ELEM_VOL REAL(DP), ALLOCATABLE :: ELEM_FIRST_MOM (:) REAL(DP), ALLOCATABLE :: ELEM_CENTROID (:) REAL(DP), ALLOCATABLE :: ELEM_SEC_MOM (:, :) REAL(DP), ALLOCATABLE :: EL_CG_SEC_MOM (:, :) !b LOGICAL, ALLOCATABLE :: ON_BOUNDARY (:) !b LOGICAL :: MASS_SYS_ALLOC = .FALSE. LOGICAL :: MASS_ELM_ALLOC = .FALSE. contains SUBROUTINE ASSEMBLE_MESH_GEOM (NODES, X) ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * ! ASSEMBLE MESH VOLUME, CENTROID, AND INERTIA MATRIX ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * Use System_Constants Use Elem_Type_Data Use Sys_Properties_Data ! for DELETED element status IMPLICIT NONE ! ALWAYS USED INTEGER, INTENT(IN) :: NODES (L_S_TOT, NOD_PER_EL) REAL(DP), INTENT(IN) :: X (MAX_NP, N_SPACE) ! Allocatable Arrays REAL(DP), ALLOCATABLE :: PT_LT (:, :), WT_LT (:) REAL(DP), ALLOCATABLE :: MASS_LT (:, :), UNIT_LT (:) REAL(DP), ALLOCATABLE :: COORD_LT (:, :) INTEGER, ALLOCATABLE :: LT_NODES (:) REAL(DP), ALLOCATABLE :: DLG_LT (:, :, :), G_LT (:, :) ! Locals !b REAL(DP), PARAMETER :: TWO_PI = 6.28318530717958623d0 REAL(DP) :: XYZ (N_SPACE) ! RADIUS IS XYZ(1) REAL(DP) :: DET ! DETERMINANT INTEGER :: IE, IP, LT, N_GP ! AXISYMMETRIC = .FALSE., UNLESS PROBLEM XYZ(1)=R IS THE RADIUS ! COORD_LT = COORDINATES OF ELEMENT GEOMETRIC NODES ! DLG_LT = PARAMETRIC DERIVATIVES OF G_LT (AT INTEGRATION PTS) ! G_LT = PARAMETRIC INTERPOLATION FOR EL TYPE GEOMETRY ! LT_NODES = THE NOD_PER_EL INCIDENCES OF THE ELEMENT ! LT = ELEMENT TYPE NUMBER ! LT_GEOM = NUMBER OF GEOMETRIC NODES ON ELEMENT ! NODES = NODAL INCIDENCES OF ALL ELEMENTS ! N_GP = NUMBER PTS TO INTEGRATE INERTIA MATRIX ! N_SPACE = PHYSICAL SPACE DIMENSIONS OF ELEMENT ! PT_LT, WT_LT = GEOMETRIC INTEGRATION POINTS & WEIGHTS ! X = COORDINATES OF SYSTEM NODES ! XYZ = COORDINATES OF POINT, IFF AXISYMMETRIC CALL INITIALIZE_GEOM_PROP IF ( I_BUG > 0 ) WRITE (N_PRT, "(/,'BEGIN GEOMETRIC ASSEMBLY')") LT = 1 ; LAST_LT = 0 ! INITIALIZE ELEMENT TYPES DO IE = 1, N_ELEMS ! (omit any boundary segment regions) !--> GET ELEMENT TYPE NUMBER IF ( N_L_TYPE > 1) LT = L_TYPE (IE) IF ( LT /= LAST_LT ) THEN ! a new type CALL GET_ELEM_TYPE_DATA (LT) ! GET CONTROLS FOR THIS TYPE LAST_LT = LT N_GP = get_numb_inertia_pts () ! get inertia rule !--> Allocate Arrays for this type's inertia arrays IF ( ALLOCATED (WT_LT) ) THEN DEALLOCATE ( UNIT_LT, MASS_LT ) DEALLOCATE ( G_LT, DLG_LT ) DEALLOCATE ( COORD_LT ) DEALLOCATE ( LT_NODES ) DEALLOCATE ( WT_LT, PT_LT) END IF ! quadrature arrays allocated ALLOCATE ( PT_LT (LT_PARM, N_GP), WT_LT (N_GP) ) ALLOCATE ( LT_NODES (LT_GEOM) ) ALLOCATE ( COORD_LT (LT_GEOM, N_SPACE) ) ALLOCATE ( DLG_LT (LT_PARM, LT_GEOM, N_GP)) ALLOCATE ( G_LT (LT_GEOM, N_GP) ) ALLOCATE ( MASS_LT (LT_GEOM, LT_GEOM), UNIT_LT (LT_GEOM)) IF ( I_BUG > 0 ) PRINT *, 'ALLOCATED GEOM_PROP ARRAYS' UNIT_LT = 1.d0 ! set type constant !--> Get geometric quadrature for this type call GET_NUMB_QUADRATURES (PT_LT, WT_LT, N_GP) ! For efficiency, fill geom arrays at each point DO IP = 1, N_GP ! EVALUATE GEOMETRIC INTERPOLATION FUNCTIONS CALL GEOM_SHAPES (PT_LT (:, IP), G_LT (:, IP)) CALL GEOM_DERIVS (PT_LT (:, IP), DLG_LT (:, :, IP)) END DO ! for type geom integration points END IF ! same type IF ( DELETED (IE) ) CYCLE ! to skip deleted element !--> EXTRACT ELEMENT GEOMETRIC NODE NUMBERS (ALL TYPES) LT_NODES = GET_ELEM_NODES (IE, LT_GEOM, NODES) !9 !--> EXTRACT ELEMENT GEOMETRIC COORDINATES (ALL TYPES) CALL ELEM_COORD (LT_GEOM, N_SPACE, X, COORD_LT, LT_NODES) !--> Build the application dependent elem metric matrix MASS_LT = 0.d0 DO IP = 1, N_GP DET = PARM_GEOM_METRIC (DLG_LT (:, :, IP), COORD_LT) IF ( AXISYMMETRIC ) THEN ! CORRECT FOR RADIUS = XYZ (1) ! GET RADIUS (X-COORDINATE) XYZ = MATMUL (G_LT (:, IP), COORD_LT) DET = DET * XYZ (1) * TWO_PI END IF ! BODY OF REVOLUTION MASS_LT = MASS_LT + DET * WT_LT(IP) & * outer_product (G_LT (:, IP), G_LT (:, IP)) END DO ! for type geom integration points !--> GENERATE ELEMENT GEOMETRIC PROPERTIES & STORE ELEM_VOL = SUM (MASS_LT) ELEM_FIRST_MOM = MATMUL (UNIT_LT (1:LT_GEOM), & MATMUL (MASS_LT, COORD_LT)) IF ( ELEM_VOL > 0.d0 ) THEN ELEM_CENTROID = ELEM_FIRST_MOM / ELEM_VOL ELSE PRINT *, 'WARNING: ASSEMBLE_MESH_GEOM, ZERO VOLUME AT ELEMENT ', IE PRINT *, 'CHECK KEYWORD remarks VALUE IN CONTROL' N_WARN = N_WARN + 1 ! INCREMENT WARNING ELEM_CENTROID = 0.d0 !b STOP 'ERROR: ASSEMBLE_MESH_GEOM, ZERO ELEMENT VOLUME' END IF ELEM_SEC_MOM = MATMUL(TRANSPOSE(COORD_LT(:, N_SPACE:1:-1)), & MATMUL(MASS_LT, COORD_LT(:, N_SPACE:1:-1))) EL_CG_SEC_MOM = ELEM_SEC_MOM - ELEM_VOL * & !b OUTER_PRODUCT (ELEM_CENTROID(N_SPACE:1:-1), & !b ELEM_CENTROID(N_SPACE:1:-1)) !b IF ( I_BUG > 0 .OR. DEBUG_INERTIA ) CALL LIST_EL_GEOM_PROP (IE) VOLUME = VOLUME + ELEM_VOL FIRST_MOMENTS = FIRST_MOMENTS + ELEM_FIRST_MOM SECOND_MOMENTS = SECOND_MOMENTS + ELEM_SEC_MOM END DO ! for all elements !--> GET SYSTEM CENTROID AND TRANSFER INERTIA THERE IF ( VOLUME > 0.d0 ) CENTROID = FIRST_MOMENTS / VOLUME CG_SECOND_MOM = SECOND_MOMENTS - VOLUME * & OUTER_PRODUCT (CENTROID(N_SPACE:1:-1), & CENTROID(N_SPACE:1:-1)) !--> GEOMETRIC ASSEMBLY COMPLETED, LIST THEM CALL LIST_SYSTEM_GEOM_PROP DEALLOCATE ( UNIT_LT, MASS_LT ) DEALLOCATE ( G_LT ) DEALLOCATE ( DLG_LT ) DEALLOCATE ( COORD_LT ) DEALLOCATE ( LT_NODES ) DEALLOCATE ( WT_LT, PT_LT) LAST_LT = 0 ! force reset of types for application IF ( I_BUG > 0 ) PRINT *, 'DEALLOCATED GEOM_PROP ARRAYS' IF ( I_BUG > 0 ) PRINT *, 'RESETTING ELEMENT TYPE' IF ( I_BUG > 0 ) WRITE (N_PRT, "(/,'END GEOMETRIC ASSEMBLY')") END SUBROUTINE ASSEMBLE_MESH_GEOM FUNCTION GET_NUMB_INERTIA_PTS () RESULT (N_GP) ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * ! RECOVER NUMBER OF POINTS TO INTEGRATE INERTIA MATRIX ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * IMPLICIT NONE INTEGER :: N_GP ! LT_SHAP = ELEMENT SHAPE CODE ! N_GP = NUMBER OF QUADRATURE POINTS ! LT_PARM = DIMENSION OF SPACE ! SELECT CASE ( LT_SHAP ) ! For max efficiency add consideration of LT_GEOM CASE (0:1) ; N_GP = 3 ! LINE ELEMENT CASE (2) ; N_GP = 7 ! TRIANGULAR ELEMENT CASE (3) ; N_GP = 9 ! QUADRILATERAL ELEMENT CASE (4) ; N_GP = 27 ! HEXAHEDRAL ELEMENT CASE (5) ; N_GP = 5 ! TETRAHEDRAL ELEMENT CASE (6) ; N_GP = 12 ! WEDGE ELEMENT CASE (7) ! USER SUPPLIED ! CALL USER_INERTIA_PT (N_GP, LT_PARM, PT, WT) STOP 'USER OPTION MISSING, GET_NUMB_INERTIA_PTS' CASE DEFAULT STOP 'INVALID OPTION, GET_NUMB_INERTIA_PTS' END SELECT END FUNCTION GET_NUMB_INERTIA_PTS subroutine ALLOCATE_GEOM_PROP ALLOCATE (FIRST_MOMENTS (N_SPACE)) ALLOCATE (CENTROID (N_SPACE)) ALLOCATE (SECOND_MOMENTS (N_SPACE, N_SPACE)) ALLOCATE (CG_SECOND_MOM (N_SPACE, N_SPACE)) ALLOCATE (ON_BOUNDARY (N_ELEMS)) ; ON_BOUNDARY = .FALSE. !b MASS_SYS_ALLOC = .TRUE. ALLOCATE (ELEM_FIRST_MOM (N_SPACE)) ALLOCATE (ELEM_CENTROID (N_SPACE)) ALLOCATE (ELEM_SEC_MOM (N_SPACE, N_SPACE)) ALLOCATE (EL_CG_SEC_MOM (N_SPACE, N_SPACE)) !b MASS_ELM_ALLOC = .TRUE. IF ( I_BUG > 0 ) PRINT *, 'ALLOCATED SYSTEM_GEOM_PROP ARRAYS' end subroutine ALLOCATE_GEOM_PROP subroutine DEALLOCATE_GEOM_PROP DEALLOCATE (EL_CG_SEC_MOM ) !b DEALLOCATE (ELEM_SEC_MOM ) DEALLOCATE (ELEM_CENTROID ) DEALLOCATE (ELEM_FIRST_MOM ) MASS_ELM_ALLOC = .FALSE. DEALLOCATE (ON_BOUNDARY ) DEALLOCATE (CG_SECOND_MOM ) DEALLOCATE (SECOND_MOMENTS ) DEALLOCATE (CENTROID ) DEALLOCATE (FIRST_MOMENTS ) MASS_SYS_ALLOC = .FALSE. IF ( I_BUG > 0 ) PRINT *, 'DEALLOCATED SYSTEM_GEOM_PROP ARRAYS' end subroutine DEALLOCATE_GEOM_PROP subroutine initialize_geom_prop ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * ! initialize geometric properties of the system ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * IMPLICIT NONE IF ( .NOT. MASS_SYS_ALLOC ) PRINT *, 'ERROR, initialize_geom_prop' VOLUME = 0.0d0 ; FIRST_MOMENTS = 0.0d0 CENTROID = 0.d0 ; SECOND_MOMENTS = 0.0d0 CG_SECOND_MOM = 0.d0 IF ( .NOT. MASS_ELM_ALLOC ) PRINT *, 'ERROR, initialize_geom_prop' ELEM_VOL = 0.0d0 ; ELEM_FIRST_MOM = 0.0d0 ELEM_SEC_MOM = 0.0d0 ; EL_CG_SEC_MOM = 0.0d0 !b end subroutine initialize_geom_prop subroutine LIST_EL_GEOM_PROP (IE) ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * ! list geometric properties of an element ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * IMPLICIT NONE INTEGER, INTENT(IN) :: IE ! Current elem number integer :: j ! loop if ( IE == 1 ) then ; write (N_PRT, *) ' ' write (N_PRT, *) '*** ELEMENT GEOMETRIC PROPERTIES ***' end if write (N_PRT, *) 'ELEM = ', IE, ', VOLUME = ', ELEM_VOL write (N_PRT, *) 'CENTROID = ', ELEM_CENTROID write (N_PRT, *) 'CENTROIDAL INERTIA MATRIX =' !b do j = 1, N_SPACE !b write (N_PRT, *) EL_CG_SEC_MOM (j, :) !b end do !b write (N_PRT, *) 'GLOBAL INERTIA MATRIX =' !b do j = 1, N_SPACE write (N_PRT, *) ELEM_SEC_MOM (j, :) end do if ( IE == N_ELEMS ) write (N_PRT, *) ' ' end subroutine LIST_EL_GEOM_PROP subroutine list_system_geom_prop ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * ! list geometric properties of the system ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * IMPLICIT NONE integer :: j ! loop write (N_PRT, *) ' ' if ( VOLUME <= 0.d0 ) then write (N_PRT, *) 'WARNING: SYSTEM GEOMETRIC PROPERTIES ARE ZERO' N_WARN = N_WARN + 1 ! INCREMENT WARNING else write (N_PRT, *) '*** SYSTEM GEOMETRIC PROPERTIES ***' write (N_PRT, "('VOLUME = ', 1PE12.5)") VOLUME if ( N_SPACE == 1 ) THEN write (N_PRT, "('CENTROID = ', 1PE12.5)") CENTROID write (N_PRT, "('INERTIA MATRIX = ', 1PE12.5)") SECOND_MOMENTS write (N_PRT, "('INERTIA MATRIX, AT CENTROID = ', 1PE12.5)") & CG_SECOND_MOM else write (N_PRT, "('CENTROID =', 3(1X, 1PE12.5))") CENTROID write (N_PRT, "(/, 'INERTIA MATRIX =')") do j = 1, N_SPACE write (N_PRT, '(3(1X, 1PE12.5))') SECOND_MOMENTS (j, 1:N_SPACE) end do write (N_PRT, "(/,'INERTIA MATRIX, AT CENTROID =')") do j = 1, N_SPACE write (N_PRT, '(3(1X, 1PE12.5))') CG_SECOND_MOM (j, 1:N_SPACE) end do write (N_PRT, *) & '(INERTIA TENSOR HAS OPPOSITE OFF-DIAGONAL SIGNS)' write (N_PRT, *) ' ' end if ! 1_D end if end subroutine list_system_geom_prop FUNCTION PARM_GEOM_METRIC (DELTA, GEOMETRY) RESULT (FFM) ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * ! FUNDAMENTAL MAGNITUDE FROM PARAMETRIC TO GEOMETRIC SPACE ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * IMPLICIT NONE REAL(DP), INTENT(IN) :: DELTA (LT_PARM, LT_GEOM) REAL(DP), INTENT(IN) :: GEOMETRY (LT_GEOM, N_SPACE) REAL(DP) :: FFM ! Automatic arrays REAL(DP) :: METRIC (LT_PARM, LT_PARM) REAL(DP) :: PARM_GEOM (LT_PARM, N_SPACE) ! GEOMETRY = COORDINATES OF THE ELEMENT'S GEOMETRIC NODES ! DELTA = LOCAL DERIVATIVES OF THE GEOMETRIC SHAPE FUNCTIONS ! FFM = DET(A), D_PHYSICAL = FFM * D_PARAMETRIC ! LT_GEOM = NUMBER OF NODES DEFINING THE GEOMETRY ! METRIC = 1ST FUND. MAG. (METRIC MATRIX) ! PARM_GEOM = PARAMETRIC DERIVATIVES OF PHYSICAL SPACE ! ESTABLISH DIFFERENTIAL LENGTHS PARM_GEOM = MATMUL (DELTA, GEOMETRY) ! FORM METRIC MATRIX METRIC = MATMUL (PARM_GEOM, transpose(PARM_GEOM)) SELECT CASE (LT_PARM) ! size of parametric space CASE (1) ; FFM = METRIC (1, 1) CASE (2) ; FFM = METRIC (1, 1) * METRIC (2, 2) & - METRIC (1, 2) * METRIC (2, 1) CASE (3) ; FFM = METRIC(1,1)*( METRIC(2,2)*METRIC(3,3) & - METRIC(3,2)*METRIC(2,3)) & + METRIC(1,2)*(-METRIC(2,1)*METRIC(3,3) & + METRIC(3,1)*METRIC(2,3)) & + METRIC(1,3)*( METRIC(2,1)*METRIC(3,2) & - METRIC(3,1)*METRIC(2,2)) CASE DEFAULT ; STOP 'INVALID LT_PARM, PARM_GEOM_METRIC' END SELECT ! LT_PARM FFM = SQRT (FFM) END FUNCTION PARM_GEOM_METRIC END MODULE GEOMETRIC_PROPERTIES Module SCP_Type_Data ! Super Convergence Patch Use System_Constants Use Elem_Type_Data IMPLICIT NONE INTEGER :: LAST_SCP, SCP_FREE INTEGER :: SCP_N, SCP_LT, SCP_GEOM, SCP_PARM, SCP_SHAP INTEGER :: REC_COUNT, SCP_RECORD_SIZE !b Integer :: SCP_FACES, SCP_SONS, SCP_SHAP, SCP_USER LOGICAL :: PATCH_ALLOC_STATUS = .FALSE. LOGICAL :: SCP_APLY_ALLOC = .FALSE. !b LOGICAL :: SCP_DATA_ALLOC = .FALSE. !b LOGICAL :: SCP_DOF_ALLOC = .FALSE. !b LOGICAL :: SCP_EQS_ALLOC = .FALSE. LOGICAL :: SCP_GAUS_ALLOC = .FALSE. !b LOGICAL :: SCP_GEOM_ALLOC = .FALSE. !b LOGICAL :: SCP_NTRP_ALLOC = .FALSE. LOGICAL :: SCP_SCAL_ALLOC = .FALSE. !b LOGICAL :: SCP_TOPO_ALLOC = .FALSE. REAL(DP), ALLOCATABLE :: SCP_H (:) ! scalar at pt REAL(DP), ALLOCATABLE :: SCP_DLH (:, :) ! parm deriv at pt REAL(DP), ALLOCATABLE :: SCP_PT (:) ! pt scp_parm coord !b ! Basic type geometry arrays !b INTEGER, ALLOCATABLE :: SCP_NODES (:) ! node numbers contains ! ---------------- BEGIN SCP RECOVERY SECTION ----------------------- SUBROUTINE LIST_SCP_ALLOC_STATUS PRINT *, ' ' PRINT *, '** SCP ALLOCATION STATUS **' PRINT *, 'SCP_APLY_ALLOC = ', SCP_APLY_ALLOC !b PRINT *, 'SCP_DATA_ALLOC = ', SCP_DATA_ALLOC !b PRINT *, 'SCP_DOF_ALLOC = ', SCP_DOF_ALLOC !b PRINT *, 'SCP_EQS_ALLOC = ', SCP_EQS_ALLOC PRINT *, 'SCP_GAUS_ALLOC = ', SCP_GAUS_ALLOC !b PRINT *, 'SCP_GEOM_ALLOC = ', SCP_GEOM_ALLOC !b PRINT *, 'SCP_NTRP_ALLOC = ', SCP_NTRP_ALLOC PRINT *, 'SCP_SCAL_ALLOC = ', SCP_SCAL_ALLOC !b PRINT *, 'SCP_TOPO_ALLOC = ', SCP_TOPO_ALLOC END SUBROUTINE LIST_SCP_ALLOC_STATUS !b SUBROUTINE ALLOCATE_SCP_APPLICATION !b ALLOCATE ( SCP_NODES (SCP_N) ) ! TOPOLOGY !b SCP_TOPO_ALLOC = .TRUE. !b IF ( I_BUG > 0 ) PRINT *, 'ALLOCATED SCP TOPOLOGY' !b END SUBROUTINE ALLOCATE_SCP_APPLICATION !b SUBROUTINE DEALLOCATE_SCP_APPLICATION !b DEALLOCATE ( SCP_NODES ) !b SCP_TOPO_ALLOC = .FALSE. !b IF ( I_BUG > 0 ) PRINT *, 'DEALLOCATED SCP TOPOLOGY' !b SCP_APLY_ALLOC = .FALSE. !b END SUBROUTINE DEALLOCATE_SCP_APPLICATION SUBROUTINE ALLOCATE_SCP_QUADRATURE IF ( .NOT. SCP_GAUS_ALLOC ) THEN ALLOCATE ( WT (LT_QP) ) ALLOCATE ( PT (LT_PARM, LT_QP) ) SCP_GAUS_ALLOC = .TRUE. IF ( I_BUG > 0 ) PRINT *, 'ALLOCATED TYPE QUADRATURE' ELSE PRINT *, 'WARNING: TYPE QP ALREADY ALLOCATED' N_WARN = N_WARN + 1 ! INCREMENT WARNING END IF END SUBROUTINE ALLOCATE_SCP_QUADRATURE SUBROUTINE DEALLOCATE_SCP_QUADRATURE IF ( SCP_GAUS_ALLOC ) THEN DEALLOCATE ( PT ) DEALLOCATE ( WT ) SCP_GAUS_ALLOC = .FALSE. IF ( I_BUG > 0 ) PRINT *, 'DEALLOCATED TYPE QUADRATURE' ELSE PRINT *, 'WARNING: TYPE QP ALREADY DEALLOCATED' N_WARN = N_WARN + 1 ! INCREMENT WARNING END IF END SUBROUTINE DEALLOCATE_SCP_QUADRATURE SUBROUTINE ALLOCATE_SCP_INTERPOLATIONS IF ( .NOT. SCP_SCAL_ALLOC ) THEN ALLOCATE ( SCP_H (SCP_N) ) ! SCALARS ALLOCATE ( SCP_DLH (SCP_PARM, SCP_N) ) ! SCALAR PARM DERIV SCP_SCAL_ALLOC = .TRUE. IF ( I_BUG > 0 ) PRINT *, 'ALLOCATED SCP SCALARS' ELSE PRINT *, 'WARNING: SCP SCAL ALREADY ALLOCATED' N_WARN = N_WARN + 1 ! INCREMENT WARNING END IF END SUBROUTINE ALLOCATE_SCP_INTERPOLATIONS SUBROUTINE DEALLOCATE_SCP_INTERPOLATIONS IF ( SCP_SCAL_ALLOC ) THEN DEALLOCATE ( SCP_DLH ) DEALLOCATE ( SCP_H ) SCP_SCAL_ALLOC = .FALSE. IF ( I_BUG > 0 ) PRINT *, 'DEALLOCATED SCP SCALARS' ELSE PRINT *, 'WARNING: SCP SCAL ALREADY DEALLOCATED' N_WARN = N_WARN + 1 ! INCREMENT WARNING END IF END SUBROUTINE DEALLOCATE_SCP_INTERPOLATIONS FUNCTION GET_SCP_ELEM_FIT (SCP_AVERAGES) RESULT (FIT_LT) ! * * * * * * * * * * * * * * * * * * * * * * * * * * ! EXTRACT ELEMENT SCP FIT FROM SYSTEM SCP FIT ! * * * * * * * * * * * * * * * * * * * * * * * * * * Use System_Constants ! for MAX_NP, SCP_FIT Use Elem_Type_Data ! for LT_N IMPLICIT NONE REAL(DP), INTENT (IN) :: SCP_AVERAGES (MAX_NP, SCP_FIT) REAL(DP) :: FIT_LT (LT_N, SCP_FIT) INTEGER :: I ! ELEM_NODES = THE LT_N NODAL INCIDENCES OF THE ELEMENT ! FIT_LT = ELEMENT NODAL AVERAGES OF GRADIENTS ! LT_N = NUMBER OF NODES FOR THIS ELEMENT TYPE ! MAX_NP = NUMBER OF NODES IN MESH ! SCP_FIT = NUMBER OF GRADIENT COMPONENTS AVERAGED ! SCP_AVERAGES = GLOBAL NODAL AVERAGES OF GRADIENTS DO I = 1, LT_N ! LOOP OVER ELEMENT NODES IF ( ELEM_NODES (I) > 0 ) THEN ! A ACTIVE NODE FIT_LT (I, :) = SCP_AVERAGES (ELEM_NODES (I), :) ELSE FIT_LT (I, :) = 0.d0 END IF END DO ! OVER NODES END FUNCTION GET_SCP_ELEM_FIT SUBROUTINE DETERMINE_SCP_BOUNDS (L_IN_PATCH, MEMBERS, NODES, X, & XYZ_MIN, XYZ_MAX, SCP_QP_SUM) ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * ! GET BOUNDS ON SCP PATCH COORDS, AND QUADRATURE POINT COUNT ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * Use System_Constants Use Elem_Type_Data ! for: ! LT_FREE, LT_GEOM, LT_N, LT_PARM, LT_QP, ELEM_NODES (LT_N), & ! COORD (LT_N, N_SPACE), GEOMETRY (LT_GEOM, N_SPACE), & Use Sys_Properties_Data Use Interface_Header IMPLICIT NONE INTEGER, INTENT (IN) :: L_IN_PATCH INTEGER, INTENT (IN) :: MEMBERS (L_IN_PATCH) INTEGER, INTENT (IN) :: NODES (L_S_TOT, NOD_PER_EL) REAL(DP), INTENT (IN) :: X (MAX_NP, N_SPACE) INTEGER, INTENT (OUT) :: SCP_QP_SUM REAL(DP), INTENT (OUT) :: XYZ_MIN (N_SPACE) REAL(DP), INTENT (OUT) :: XYZ_MAX (N_SPACE) REAL(DP) :: EL_MIN (N_SPACE), EL_MAX (N_SPACE) REAL(DP) :: BORDER (N_SPACE) INTEGER :: IE, IE_FIRST, IM, LT ! VARIABLES: ! DELETED = TRUE IF THE ELEMENT HAS BEEN DELETED ! ELEM_NODES = THE NOD_PER_EL INCIDENCES OF THE ELEMENT ! INDEX = SYSTEM DOF NUMBERS ASSOCIATED WITH ELEMENT ! L_HOMO = 1, IF ELEMENT PROPERTIES ARE HOMOGENEOUS ! L_S_TOT = TOTAL NUMBER OF ELEMENTS & THEIR SEGMENTS ! L_SHAPE = SHAPE FLAG FOR QUADRATURE RULE SELECTION ! LT_DATA = ARRAY OF CONSTANTS FOR EACH ELEMENT TYPE ! LT = ELEMENT TYPE NUMBER (IF USED) ! LT_QP = NUMBER OF QUADRATURE PTS FOR ELEMENT TYPE ! LT_USER = OPTIONAL USER CODE ASSOCIATED WITH TYPE ! L_TYPE = ELEMENT TYPE NUMBER ! MAX_NP = NUMBER OF SYSTEM NODES ! NOD_PER_EL = NUMBER OF NODES PER ELEMENT ! NODES = NODAL INCIDENCES OF ALL ELEMENTS ! N_QP = NUMBER OF QUADRATURE POINTS, >= LT_QP ! N_SPACE = DIMENSION OF SPACE ! X = COORDINATES OF SYSTEM NODES ! XYZ = SPACE COORDINATES AT A POINT XYZ_MIN = HUGE (1.d0) ; XYZ_MAX = -HUGE (1.d0) ! INITIALIZE EL_MIN = HUGE (1.d0) ; EL_MAX = -HUGE (1.d0) ! INITIALIZE SCP_QP_SUM = 0 ; LT = 1 ! INITIALIZE ! THE FIRST MEMBER DEFINES THE PATCH DETAILS IE_FIRST = MEMBERS (1) IF ( IE_FIRST == 0 ) RETURN ! INACTIVE ELEMENT IF ( DELETED (IE_FIRST) ) RETURN ! INACTIVE ELEMENT !--> GET FIRST ELEMENT TYPE NUMBER IF ( N_L_TYPE > 1) LT = L_TYPE (IE_FIRST) ! SAME AS LAST TYPE ? IF ( LT /= LAST_LT ) THEN ! this is a new type CALL GET_ELEM_TYPE_DATA (LT) ! GET CONTROLS FOR ELEM TYPE LAST_LT = LT !b following moved to CALC_SCP_AVE_NODE_FLUXES may 2001 !b CALL SET_SCP_TYPE_DATA (LT) ! GET CONTROLS FOR SCP TYPE CALL SET_SCP_TYPE_DATA (LT) ! GET CONTROLS FOR SCP TYPE END IF ! a new element type and scp type !---- DO IM = 1, L_IN_PATCH IE = MEMBERS (IM) IF ( IE /= 0 ) THEN ! ACTIVE ELEMENT !--> GET ELEMENT TYPE NUMBER IF ( N_L_TYPE > 1) LT = L_TYPE (IE) ! SAME AS LAST TYPE ? !b IF ( LT /= LAST_LT ) CALL GET_ELEM_TYPE_DATA (LT) IF ( LT /= LAST_LT ) THEN CALL GET_ELEM_TYPE_DATA (LT) CALL DEALLOCATE_TYPE_APPLICATION CALL ALLOCATE_TYPE_APPLICATION LAST_LT = LT IF ( DEBUG_EL_TYPE ) THEN CALL LIST_TYPE_ALLOC_STATUS PRINT *,'IE, SIZE (ELEM_NODES) ', IE, SIZE (ELEM_NODES) END IF ! debug END IF ! TYPE CHANGE IF ( DELETED (IE) ) CYCLE ! to skip deleted element SCP_QP_SUM = SCP_QP_SUM + LT_QP ! UPDATE QUADRATURE COUNT ! EXTRACT ELEMENT NODE NUMBERS & COORDINATES IF ( DEBUG_SCP ) THEN PRINT *,' IM, IE, LT, LT_N ', IM, IE, LT, LT_N PRINT *,'SIZE(ELEM_NODES) ', SIZE(ELEM_NODES) END IF ! XXX ? problem with patch of mixed element types ? !b ELEM_NODES = GET_ELEM_NODES (IE, LT_N, NODES) ELEM_NODES (1:LT_N) = GET_ELEM_NODES (IE, LT_N, NODES) !b CALL ELEM_COORD (LT_N, N_SPACE, X, COORD, ELEM_NODES) !b EL_MAX (:) = MAXVAL ( COORD (1:LT_N, :) ) ! bug fix !b EL_MIN (:) = MINVAL ( COORD (1:LT_N, :) ) ! bug fix EL_MAX (:) = MAXVAL ( COORD (1:LT_N, :), DIM = 1 ) ! bug fix EL_MIN (:) = MINVAL ( COORD (1:LT_N, :), DIM = 1 ) ! bug fix WHERE ( EL_MAX > XYZ_MAX ) XYZ_MAX = EL_MAX WHERE ( EL_MIN < XYZ_MIN ) XYZ_MIN = EL_MIN END IF ! ACTIVE ELEMENT END DO ! OVER IM PATCH MEMBERS ! VALIDATE PATCH SHAPE IF ( ANY( (XYZ_MAX - XYZ_MIN) == 0.D0 )) THEN PRINT *,'ERROR, PATCH ', IE_FIRST, ' HAS A ZERO DIMENSION' PRINT *,'MAX BOUNDS ARE: ', XYZ_MAX PRINT *,'MIN BOUNDS ARE: ', XYZ_MIN STOP 'INVALID PATCH IN DETERMINE_SCP_BOUNDS' END IF ! FINITE SHAPED PATCH ! ADD A SMALL BORDER TO PATCH BORDER = ABS ( XYZ_MAX - XYZ_MIN ) * 0.01 XYZ_MAX = XYZ_MAX + BORDER XYZ_MIN = XYZ_MIN - BORDER END SUBROUTINE DETERMINE_SCP_BOUNDS FUNCTION GET_SCP_PT_AT_XYZ (XYZ, XYZ_MIN, XYZ_MAX) RESULT (POINT) ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * ! FIND LOCAL POINT IN PATCH FOR GIVEN XYZ POSITION ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * Use System_Constants Use Elem_Type_Data IMPLICIT NONE REAL(DP), INTENT (IN) :: XYZ (N_SPACE) REAL(DP), INTENT (IN) :: XYZ_MIN (N_SPACE), XYZ_MAX (N_SPACE) REAL(DP) :: POINT (N_SPACE) LOGICAL :: FATAL ! debug aid !b ! . S PATCH MAPPING: ELEMENT AND ITS NEIGHBORS ARE BOUNDED BY ! . XYZ_MIN < XYZ < XYZ_MAX, WHICH MATCHES THE NATURAL PATCH ! |\ AT -1 < ABC < 1, OR THE FULLY INTERIOR REGION OF THE ! | \ UNIT PATCH AT 0 < RST < (1/N_SPACE). THIS MEANS THAT ! | \ MAKING THE PATCH AXES PARALLEL TO XYZ AXES WILL CAUSE ! | \ A CONSTANT DIAGONAL JACOBIAN BETWEEN XYZ AND THE PATCH. ! | \ ! |----------\ XYZ_MAX, (1,1), (1/M, 1/M); M = N_SPACE ! | .B | \ ! | . | \ THUS, COMPUTING THE LOCAL COORDINATE OF A SAMPLE ! | ...A | \ POINT IS A SIMPLE LINEAR TRANSFORMATION. ! | *P | \ NOTE: 0 <= SUM(RST) <= 1. ! |----------|---------\ .....R ! XYZ_MIN, (-1,-1), (0,0) [ONLY WEDGE FAMILY CURRENTLY EXCLUDED.] IF ( SCP_PARM /= N_SPACE ) STOP & 'INVALID SCP TRANSFORMATION TYPE, GET_SCP_PT_AT_XYZ' IF (SCP_SHAP == 1 .OR. SCP_SHAP == 3 .OR. SCP_SHAP == 4 ) THEN ! A NATURAL COORDINATE PATCH TRANSFORMATION, ! FOR A: LINE, QUADRILATERAL, OR HEXAHEDRAL PATCH POINT = (2.D0 * XYZ - (XYZ_MAX + XYZ_MIN)) / (XYZ_MAX - XYZ_MIN) ELSE IF (SCP_SHAP == 2 .OR. SCP_SHAP == 5 ) THEN ! A UNIT COORDINATE PATCH TRANSFORMATION, ! FOR A: TRIANGULAR, OR TETRAHEDRAL PATCH POINT = (XYZ - XYZ_MIN) / (SCP_PARM * (XYZ_MAX - XYZ_MIN)) ELSE PRINT *, 'ERROR, NO CONVERSION FOR SHAPE ', SCP_SHAP STOP 'SHAPE ERROR IN GET_SCP_PT_AT_XYZ' END IF ! PARAMETRIC PARENT PATCH SHAPE' ! VALIDATE POINT LIES IN PARENT SPACE FATAL = .FALSE. IF (SCP_SHAP == 2 .OR. SCP_SHAP == 5 ) THEN ! A UNIT COORDINATE PATCH TRANSFORMATION IF ( ANY (POINT < 0.d0) ) FATAL = .TRUE. !b !b STOP 'PATCH MAP OUT OF BOUNDS, GET_SCP_PT_AT_XYZ' IF ( SUM (POINT) > 1.d0 ) FATAL = .TRUE. !b !b STOP 'PATCH MAP OUT OF BOUNDS, GET_SCP_PT_AT_XYZ' ELSE ! A NATURAL COORDINATE PATCH TRANSFORMATION IF ( ANY (ABS(POINT) > 1.d0) ) FATAL = .TRUE. !b !b STOP 'PATCH MAP OUT OF BOUNDS, GET_SCP_PT_AT_XYZ' END IF IF ( FATAL ) THEN PRINT *,'ERROR: PATCH MAP OUT OF BOUNDS, GET_SCP_PT_AT_XYZ' PRINT *,'FOR GLOBAL XYZ = ', XYZ PRINT *,'FOR PARAMETRIC RS = ', POINT STOP 'PATCH MAP OUT OF BOUNDS, GET_SCP_PT_AT_XYZ' END IF ! impossible geometry END FUNCTION GET_SCP_PT_AT_XYZ !SUBROUTINE SET_SCP_TYPE_DATA_OLD (LT) !! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * !! GET DATA BASED ON ELEMENT TYPE (SEE CONTROL_INPUT ) !! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * !Use System_Constants !Use Elem_Type_Data ! for LAST_LT, N_L_TYPE ! IMPLICIT NONE ! INTEGER, INTENT (IN) :: LT ! !b INTEGER :: LT_FAKE ! !! LAST_LT = LAST ACTIVE ELEMENT TYPE NUMBER !! LAST_SCP = LAST ACTIVE PATCH TYPE NUMBER !! LT = ELEMENT TYPE NUMBER !! LT_N = NUMBER OF NODES PER ELEMENT !! LT_GEOM = NUMBER OF ELEMENT GEOMERTY NODES !! LT_PARM = NUMBER OF PARAMETRIC SPACES FOR ELEMENT !! LT_SHAP = ELEMENT SHAPE FLAG NUMBER !! N_L_TYPE = NUMBER OF ELEMENT TYPES !! SCP_LT = PATCH TYPE NUMBER !! SCP_N = NUMBER OF NODES PER PATCH !! SCP_GEOM = NUMBER OF PATCH GEOMERTY NODES !! SCP_PARM = NUMBER OF PARAMETRIC SPACES FOR PATCH !! SCP_SHAP = PATCH SHAPE FLAG NUMBER ! ! IF ( LT < 1 .OR. LT > N_L_TYPE ) STOP & ! 'ELEMENT TYPE WRONG, SET_SCP_TYPE_DATA_OLD' ! !! GET CONTROLS FOR THIS PATCH TYPE ! SELECT CASE ( SCP_INC ) ! increase in degree ! CASE (0) ! same as reference element ! SCP_LT = LT ! SCP_N = LT_N ! SCP_GEOM = LT_GEOM !b un-used ?? ! SCP_PARM = LT_PARM ! SCP_SHAP = LT_SHAP ! ! 1 or 2 add logic here !b xxx ! CASE DEFAULT ! STOP 'SET_SCP_TYPE_DATA_OLD: NO DEGREE INCREASE ALLOWED' ! END SELECT !END SUBROUTINE SET_SCP_TYPE_DATA_OLD SUBROUTINE SET_SCP_TYPE_DATA (LT) ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * ! GET DATA BASED ON ELEMENT TYPE (SEE CONTROL_INPUT ) ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * Use System_Constants Use Elem_Type_Data ! for LAST_LT, N_L_TYPE IMPLICIT NONE INTEGER, INTENT (IN) :: LT !b INTEGER :: LT_FAKE ! LAST_LT = LAST ACTIVE ELEMENT TYPE NUMBER ! LAST_SCP = LAST ACTIVE PATCH TYPE NUMBER ! LT = ELEMENT TYPE NUMBER ! LT_N = NUMBER OF NODES PER ELEMENT ! LT_GEOM = NUMBER OF ELEMENT GEOMERTY NODES ! LT_PARM = NUMBER OF PARAMETRIC SPACES FOR ELEMENT ! LT_QP = NUMBER OF QUADRATURE POINTS FOR ELEMENT TYPE ! LT_SHAP = ELEMENT SHAPE FLAG NUMBER ! N_L_TYPE = NUMBER OF ELEMENT TYPES ! SCP_LT = PATCH TYPE NUMBER ! SCP_N = NUMBER OF NODES PER PATCH ! SCP_GEOM = NUMBER OF PATCH GEOMERTY NODES ! SCP_PARM = NUMBER OF PARAMETRIC SPACES FOR PATCH ! SCP_QP = NUMBER OF QUADRATURE POINTS NEEDED IN A SCP ! SCP_SHAP = PATCH SHAPE FLAG NUMBER !b ! SPECIAL CASE OF FITTING DOF FIRST !b IF ( SCP_DOF ) THEN ! use quadrilaterals only !b END IF ! special case IF ( LT < 1 .OR. LT > N_L_TYPE ) STOP & 'ELEMENT TYPE WRONG, SET_SCP_TYPE_DATA' ! GET CONTROLS FOR THIS PATCH TYPE SCP_LT = LT SCP_PARM = LT_PARM SCP_SHAP = LT_SHAP SCP_GEOM = LT_GEOM !b un-used ?? SELECT CASE ( SCP_INC ) ! increase in degree CASE (0) ! same as reference element SCP_N = LT_N CALL GET_PATCH_QUADRATURE_ORDER (LT_SHAP, LT_QP, SCP_QP) CASE (1) ! increase patch polynomial degree by 1 CALL GET_PATCH_QUADRATURE_ORDER (LT_SHAP, LT_QP, SCP_QP) SELECT CASE ( LT_SHAP ) CASE (1) ! LINE ELEMENT IF ( LT_N < 4 ) THEN SCP_N = LT_N + 1 ELSE SCP_N = LT_N PRINT *, 'WARNING, NO PATCH DEGREE INCREASE' END IF ! current lib limit CASE (2) ! TRIANGULAR ELEMENT SELECT CASE ( NOD_PER_EL ) CASE ( 3) ; SCP_N = 6 CASE ( 6) ; SCP_N = 10 CASE (10) ; SCP_N = 15 ! Current limit CASE (15) CASE DEFAULT PRINT *, 'WARNING, NO PATCH DEGREE INCREASE' SCP_N = 15 END SELECT ! number of element nodes CASE (3) ! QUADRILATERAL ELEMENT SELECT CASE ( NOD_PER_EL ) CASE ( 4) ; SCP_N = 9 CASE ( 8) ; SCP_N = 9 CASE ( 9) ; SCP_N = 17 ! Current limit CASE (17) CASE DEFAULT PRINT *, 'WARNING, NO PATCH DEGREE INCREASE' SCP_N = 17 END SELECT ! number of element nodes CASE (4) ! HEXAHEDRAL ELEMENT SELECT CASE ( NOD_PER_EL ) CASE ( 8) ; SCP_N = 20 ! Current limit CASE ( 8) CASE DEFAULT PRINT *, 'WARNING, NO PATCH DEGREE INCREASE' SCP_N = 20 END SELECT ! number of element nodes CASE (5) ! TETRAHEDRAL ELEMENT ! Current limit CASE ( 4) PRINT *, 'WARNING, NO PATCH DEGREE INCREASE' SCP_N = 4 !b CASE (6) ! WEDGE ELEMENT CASE (7) ! USER SUPPLIED PRINT *, 'WARNING, NO PATCH DEGREE INCREASE' SCP_N = LT_N CASE DEFAULT STOP 'INVALID OPTION, SET_SCP_TYPE_DATA' END SELECT CASE (2) ! increase patch polynomial degree by 2 CALL GET_PATCH_QUADRATURE_ORDER (LT_SHAP, LT_QP, SCP_QP) CALL GET_PATCH_QUADRATURE_ORDER (LT_SHAP, SCP_QP, SCP_QP) SELECT CASE ( LT_SHAP ) CASE (1) ! LINE ELEMENT IF ( LT_N < 3 ) THEN SCP_N = LT_N + 2 ELSE SCP_N = LT_N PRINT *, 'WARNING, NO PATCH DEGREE INCREASE' END IF ! current lib limit CASE (2) ! TRIANGULAR ELEMENT SELECT CASE ( NOD_PER_EL ) CASE ( 3) ; SCP_N = 10 CASE ( 6) ; SCP_N = 15 CASE (10) ; SCP_N = 15 PRINT *, 'WARNING, PATCH DEGREE INCREASE = 1' ! Current limit CASE (15) CASE DEFAULT PRINT *, 'WARNING, NO PATCH DEGREE INCREASE' SCP_N = 15 END SELECT ! number of element nodes CASE (3) ! QUADRILATERAL ELEMENT SELECT CASE ( NOD_PER_EL ) CASE ( 4) ; SCP_N = 9 PRINT *, 'WARNING, PATCH DEGREE INCREASE = 1' CASE ( 8) ; SCP_N = 17 CASE ( 9) ; SCP_N = 17 ! Current limit CASE (17) CASE DEFAULT PRINT *, 'WARNING, PATCH DEGREE INCREASE = 1' SCP_N = 17 END SELECT ! number of element nodes CASE (4) ! HEXAHEDRAL ELEMENT SELECT CASE ( NOD_PER_EL ) CASE ( 8) ; SCP_N = 20 ! Current limit CASE ( 8) CASE DEFAULT PRINT *, 'WARNING, PATCH DEGREE INCREASE = 1' SCP_N = 20 END SELECT ! number of element nodes CASE (5) ! TETRAHEDRAL ELEMENT ! Current limit CASE ( 4) PRINT *, 'WARNING, NO PATCH DEGREE INCREASE' SCP_N = 4 !b CASE (6) ! WEDGE ELEMENT CASE (7) ! USER SUPPLIED PRINT *, 'WARNING, NO PATCH DEGREE INCREASE' SCP_N = LT_N CASE DEFAULT STOP 'INVALID OPTION, SET_SCP_TYPE_DATA' END SELECT CASE DEFAULT PRINT *,'WARNING, SET_SCP_TYPE_DATA: SET scp_deg_inc ', & 'TO 1 OR 2 IN KEYWORDS' PRINT *, 'WARNING, NO PATCH DEGREE INCREASE' N_WARN = N_WARN + 1 END SELECT END SUBROUTINE SET_SCP_TYPE_DATA ! FUNCTION GET_PATCH_DEGREE (L_IN_PATCH, MEMBERS) RESULT (DEGREE) ! Use System_Constants ! for DEG_HOMO, SCP_DEG ! Use Elem_Type_Data ! INTEGER, INTENT (IN) :: L_IN_PATCH ! NUMBER OF ELEMENTS IN PATCH ! INTEGER, INTENT (IN) :: MEMBERS (L_IN_PATCH) ! ELEMENTS IN PATCH ! INTEGER :: IM, MEMB_DEG ! LOOPS ! INTEGER :: DEGREE ! IF ( DEG_HOMO == 1 ) THEN ! HOMOGENEOUS DEGREE IN EACH PATCH ! DEGREE = SCP_DEG ! ELSE ! SEARCH FOR ELEMENT OF HIGHEST DEGREE IN MEMBERS ! DEGREE = 1 ! DEFAULT ! DO IM = 1, L_IN_PATCH ! LOOP OVER ELEMENTS ! MEMB_DEG = GET_ELEMENT_DEGREE ( MEMBERS (IM) ) ! IF ( MEMB_DEG > DEGREE ) DEGREE = MEMB_DEG ! END DO ! OVER PATCH ELEMENT MEMBERS ! END IF ! HOMOGENEOUS ELEMENT DEGREES ! END FUNCTION GET_PATCH_DEGREE ! FUNCTION GET_SCP_DEG () RESULT (DEGREE) ! Use System_Constants ! INTEGER :: DEGREE ! DEGREE = SCP_DEG ! END FUNCTION GET_SCP_DEG ! SUBROUTINE SET_SCP_DEG (DEGREE) ! Use System_Constants ! INTEGER, INTENT (IN) :: DEGREE ! SCP_DEG = DEGREE ! END SUBROUTINE SET_SCP_DEG ! ---------------- END SCP RECOVERY SECTION ----------------------- ! ----------- SINGULAR VALUE DECOMPOSITION SECTION ---------------- SUBROUTINE SVDC_FACTOR (A, M, N, W, V) ! * * * * * * * * * * * * * * * * * * * * * * * * * * ! SINGULAR VALUE DECOMPOSITION FACTORIZATION: ! GIVEN M BY N MATRIX A, DECOMPOSE IT INTO THE PRODUCTS ! A = U * W * V_TRANSPOSE, WHERE U OVERWRITES A AND IS ! RETURNED IN ITS SPACE, W IS A DIAGONAL MATRIX RETURNED ! AS A VECTOR, AND V IS A SQUARE N BY N MATRIX (WHICH IS ! RETURNED INSTEAD OF THE TRANSPOSE OF V) ! * * * * * * * * * * * * * * * * * * * * * * * * * * ! SEE NUMERICAL RECIPES IN FORTRAN 90, PRESS ET AL USE PRECISION_MODULE IMPLICIT NONE INTEGER, INTENT (IN) :: M, N REAL (DP), INTENT (INOUT) :: A (M, N) REAL (DP), INTENT (OUT) :: W (N) REAL (DP), INTENT (OUT) :: V (N, N) REAL (DP) :: RV1 (N) ! Automatic array REAL (DP) :: G, SCALE, ANORM, S, F, H, C, X, Y, Z INTEGER :: ITS, I, J, JJ, K, L, NM INTEGER, PARAMETER :: MAX_ITS = 30 IF (M < N) THEN WRITE (6, *) 'FATAL ERROR', M, N WRITE (6, *) 'Fewer equations than unknowns, SVDC_FACTOR' STOP 'Fewer equations than unknowns, SVDC_FACTOR' END IF ! Householder reduction to bidiagonal form G = 0.d0 SCALE = 0.d0 ANORM = 0.d0 DO I = 1,N L = I+1 RV1(I) = SCALE*G G = 0.d0 ! S = 0.d0 SCALE = 0.d0 IF (I <= M) THEN ! DO K = I,M ! SCALE = SCALE+ABS(A(K,I)) ! END DO SCALE = SUM ( ABS(A(I:M, I))) IF (SCALE /= 0.d0) THEN !9 use epsilon ! DO K = I,M ! A(K,I) = A(K,I)/SCALE ! S = S+A(K,I)*A(K,I) ! END DO A(I:M, I) = A(I:M, I) / SCALE S = DOT_PRODUCT (A(I:M, I), A(I:M, I)) F = A(I,I) G = -SIGN(SQRT(S),F) H = F*G-S A(I,I) = F-G IF (I /= N) THEN DO J = L,N S = 0.d0 DO K = I,M S = S+A(K,I)*A(K,J) END DO F = S/H DO K = I,M A(K,J) = A(K,J)+F*A(K,I) END DO END DO END IF DO K = I,M A(K,I) = SCALE*A(K,I) END DO END IF END IF W(I) = SCALE*G G = 0.d0 S = 0.d0 SCALE = 0.d0 IF ((I <= M) .AND. (I /= N)) THEN ! DO K = L,N ! SCALE = SCALE+ABS(A(I,K)) ! END DO SCALE = SUM ( ABS (A(I, L:N))) IF (SCALE /= 0.d0) THEN DO K = L,N A(I,K) = A(I,K)/SCALE S = S+A(I,K)*A(I,K) END DO F = A(I,L) G = -SIGN(SQRT(S),F) H = F*G-S A(I,L) = F-G DO K = L,N RV1(K) = A(I,K)/H END DO IF (I /= M) THEN DO J = L,M S = 0.d0 DO K = L,N S = S+A(J,K)*A(I,K) END DO DO K = L,N A(J,K) = A(J,K)+S*RV1(K) END DO END DO END IF DO K = L,N A(I,K) = SCALE*A(I,K) END DO END IF END IF ANORM = MAX(ANORM,(ABS(W(I))+ABS(RV1(I)))) END DO ! Accumulation of right-hand transformations DO I = N,1,-1 IF (I < N) THEN IF (G /= 0.d0) THEN DO J = L,N V(J,I) = (A(I,J)/A(I,L))/G END DO DO J = L,N S = 0.d0 DO K = L,N S = S+A(I,K)*V(K,J) END DO DO K = L,N V(K,J) = V(K,J)+S*V(K,I) END DO END DO END IF DO J = L,N V(I,J) = 0.d0 V(J,I) = 0.d0 END DO END IF V(I,I) = 1.d0 G = RV1(I) L = I END DO ! Accumulation of left-hand transformations DO I = N,1,-1 L = I+1 G = W(I) IF (I < N) THEN DO J = L,N A(I,J) = 0.d0 END DO END IF IF (G /= 0.d0) THEN G = 1.d0/G IF (I /= N) THEN DO J = L,N S = 0.d0 DO K = L,M S = S+A(K,I)*A(K,J) END DO F = (S/A(I,I))*G DO K = I,M A(K,J) = A(K,J)+F*A(K,I) END DO END DO END IF DO J = I,M A(J,I) = A(J,I)*G END DO ELSE DO J = I,M A(J,I) = 0.d0 END DO END IF A(I,I) = A(I,I)+1.d0 END DO ! Diagonalization of the bidiagonal form DO K = N,1,-1 DO ITS = 1, MAX_ITS DO L = K,1,-1 NM = L-1 IF ((ABS(RV1(L))+ANORM) == ANORM) GOTO 2 IF ((ABS(W(NM))+ANORM) == ANORM) GOTO 1 END DO 1 CONTINUE C = 0.d0 S = 1.d0 DO I = L,K F = S*RV1(I) RV1(I) = C*RV1(I) IF ((ABS(F)+ANORM) == ANORM) GOTO 2 G = W(I) H = SQRT(F*F+G*G) W(I) = H H = 1.d0/H C = G*H S = -(F*H) DO J = 1,M Y = A(J,NM) Z = A(J,I) A(J,NM) = (Y*C)+(Z*S) A(J,I) = -(Y*S)+(Z*C) END DO END DO 2 CONTINUE Z = W(K) IF (L == K) THEN IF (Z < 0.d0) THEN W(K) = -Z DO J = 1,N V(J,K) = -V(J,K) END DO END IF GOTO 3 END IF IF (ITS == MAX_ITS) WRITE (6, *) & 'No convergence in ', MAX_ITS, ' iterations' X = W(L) NM = K-1 Y = W(NM) G = RV1(NM) H = RV1(K) F = ((Y-Z)*(Y+Z)+(G-H)*(G+H))/(2.d0*H*Y) G = SQRT(F*F+1.d0) F = ((X-Z)*(X+Z)+H*((Y/(F+SIGN(G,F)))-H))/X ! Next QR transformation C = 1.d0 S = 1.d0 DO J = L,NM I = J+1 G = RV1(I) Y = W(I) H = S*G G = C*G Z = SQRT(F*F+H*H) RV1(J) = Z C = F/Z S = H/Z F = (X*C)+(G*S) G = -(X*S)+(G*C) H = Y*S Y = Y*C DO JJ = 1,N X = V(JJ,J) Z = V(JJ,I) V(JJ,J) = (X*C)+(Z*S) V(JJ,I) = -(X*S)+(Z*C) END DO Z = SQRT(F*F+H*H) W(J) = Z IF (Z /= 0.d0) THEN Z = 1.d0/Z C = F*Z S = H*Z END IF F = (C*G)+(S*Y) X = -(S*G)+(C*Y) DO JJ = 1,M Y = A(JJ,J) Z = A(JJ,I) A(JJ,J) = (Y*C)+(Z*S) A(JJ,I) = -(Y*S)+(Z*C) END DO END DO RV1(L) = 0.d0 RV1(K) = F W(K) = X END DO 3 CONTINUE END DO END SUBROUTINE SVDC_FACTOR SUBROUTINE SVDC_BACK_SUBST (U, W, V, M, N, B, X) ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * ! SOLVES A * X = B FOR VECTOR X, WHERE RECTANGULAR MATRIX A HAS ! BEEN DECOMPOSED INTO U, W, V BY SVDC_FACTOR. MAY BE CALLED ! SEQUENTIALLY WITH DIFFERENT B VALUES FOR NEW X VALUES ! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * USE Precision_Module IMPLICIT NONE INTEGER, INTENT (IN) :: M, N REAL (DP), INTENT (IN) :: U (M, N) REAL (DP), INTENT (IN) :: W (N) REAL (DP), INTENT (IN) :: V (N, N) REAL (DP), INTENT (IN) :: B (M) REAL (DP), INTENT (OUT) :: X (N) REAL (DP) :: TMP (N) ! Automatic work space WHERE ( W /= 0.D0 ) TMP = MATMUL (B, U)/W ! DIAG(1/W_J)*U_T*B ELSEWHERE TMP = 0.D0 END WHERE X = MATMUL (V, TMP) END SUBROUTINE SVDC_BACK_SUBST ! -------------- END SINGULAR VALUE DECOMPOSITION SECTION ------------- END module SCP_Type_Data MODULE keyword_buffer ! For nonadvancing keword, data input ! Copyright J. E. Akin, 1998 ! Input lines (of up to LINE_SIZE characters) have a keyword (of up ! to MAX_KEY characters) and may be followed by a mixture of ! integers, reals, or strings separated by a gap. ! Strings with more than one word should be in quotes ("), or use ! the underscore (_) instead of gaps. ! A gap is a blank space, comma, or horizontal tab. implicit none public :: get_int, get_key, get_real, get_string ! subroutines public :: MAX_KEY, LINE_SIZE, WORD_SIZE ! size parameters public :: ECHO_UNIT, ECHO_KEY, ECHO_OPEN ! for debug private ! everything not listed as public integer, parameter :: LINE_SIZE = 80, MAX_KEY = 16 integer, parameter :: WORD_SIZE = LINE_SIZE - MAX_KEY integer, parameter :: ECHO_UNIT = 17 character(len=5), parameter :: FMT_LINE = "(a80)" character(len=1), parameter :: BLA = ' ', COM = ',' character(len=1), parameter :: TAB = achar(9) character(len=LINE_SIZE) :: buffer ! input buffer line logical :: ECHO_KEY = .FALSE. ! true for debug logical :: ECHO_OPEN = .FALSE. ! echo file status integer :: LINES_READ = 0 ! line information logical :: LINE_FINISHED = .FALSE. ! line information integer :: ioResult ! 0=no error, <0=EOF integer :: NEXT_START, START, STOP, LENGTH, N_CHARS contains ! encapsulated functionality subroutine open_echo_unit !--------------------------------------------------------- ! Open the file to echo keywords that are parsed !--------------------------------------------------------- implicit none integer :: ok ! file status !b if ( ECHO_KEY ) then open (ECHO_UNIT, file = 'keyword_input.echo', & action ='write', status = 'unknown', & iostat = ok) if ( ok == 0 ) ECHO_OPEN = .TRUE. !b else !b print *,'Note: set keyword_buffer ECHO_KEY = .true.' !b ECHO_OPEN = .FALSE. !b end if end subroutine open_echo_unit function is_a_gap (new) result (gap) !--------------------------------------------------------- ! Test if a character is a blank, comma, or tab !--------------------------------------------------------- implicit none character(len=1), intent(in) :: new logical :: gap gap = .false. ! initialize if (new == BLA .or. new == TAB .or. new == COM) gap = .true. end function is_a_gap subroutine locate_next_word !--------------------------------------------------------- ! Locate start and stop characters of the next word !--------------------------------------------------------- implicit none integer :: j ! loops do j = NEXT_START, LENGTH ! find non gap if ( .not. is_a_gap (buffer(j:j)) ) then ! word starts START = j ; exit ! this loop end if end do ! for START STOP = START ! initialize NEXT_START = 1 ! initialize LINE_FINISHED = .TRUE. ! initialize do j = START, LENGTH ! find next gap if ( is_a_gap (buffer(j:j)) ) then ! found gap if ( j /= LENGTH ) LINE_FINISHED = .FALSE. NEXT_START = j ! passed word exit ! this loop else ; STOP = j ! update stop end if end do ! for STOP N_CHARS = STOP - START + 1 ! number of characters end subroutine locate_next_word subroutine get_key (key_word) !--------------------------------------------------------- ! get the first keyword on line !--------------------------------------------------------- implicit none character(len=MAX_KEY), intent(out) :: key_word character(len=13) :: S_FORMAT ! for debug ! read a line of input into the buffer key_word = ' ' read (5, FMT_LINE, iostat = ioResult) buffer ! input buffer line if ( ioResult < 0 ) then ! test for EOF key_word = 'EOF' print *, 'get_key: EOF found after line ', LINES_READ ! debug else LINES_READ = LINES_READ + 1 ! increment count LENGTH = MAX (1, len_trim (buffer)) ! true length START = 1 ; STOP = LENGTH ; NEXT_START = 1 ! initialize word call locate_next_word ! find keyword key_word (1:N_CHARS) = buffer(START:STOP) ! copy keyword !K if ( ECHO_KEY ) then ! begin debug echo ! create format via internal file write (S_FORMAT, '("(A", I2, ",1X)")' ) N_CHARS write (ECHO_UNIT, S_FORMAT, advance='no') & key_word (1:N_CHARS) ! append key_word !K end if ! echo !K write (6, S_FORMAT, advance='no') & !b !K key_word (1:N_CHARS) ! append key_word !b print *, buffer (1:LENGTH) !K !b end if ! EOF end subroutine get_key subroutine get_string (string) !--------------------------------------------------------- ! get the next string on the line !--------------------------------------------------------- implicit none character(len=WORD_SIZE), intent(out) :: string character(len=13) :: S_FORMAT ! for debug integer :: j ! loops string = ' ' ; N_CHARS = 1 ! initialize if ( .not. LINE_FINISHED ) then ! a string should be present call locate_next_word ! find the string size if ( N_CHARS > 1 .and. buffer(START:START) == '"' ) then do j = START + 1, LENGTH ! locate end of quotes STOP = j if ( buffer(j:j) == '"' ) exit ! this loop end do ! for next quote N_CHARS = STOP - START + 1 ! revise number of characters end if ! string starts with a quote string (1:N_CHARS) = buffer(START:STOP) ! copy from buffer !K if ( ECHO_KEY ) then ! begin debug echo ! create format via internal file write (S_FORMAT, '("(A", I2, ",1X)")' ) N_CHARS write (ECHO_UNIT, S_FORMAT, advance='no') & string (1:N_CHARS) ! append string !K end if ! echo end if ! LINE_FINISHED end subroutine get_string subroutine get_int (fixed) !--------------------------------------------------------- ! get the next integer !--------------------------------------------------------- implicit none integer, intent(out) :: fixed ! the integer if ( LINE_FINISHED ) then ! assign default of zero fixed = 0 ; N_CHARS = 1 else ! should be present call locate_next_word ! locate integer read (buffer(START:STOP), *) fixed ! read integer from buffer !K if ( ECHO_KEY ) write (ECHO_UNIT, '(I9)', advance='no') fixed write (ECHO_UNIT, '(I9)', advance='no') fixed !K end if ! LINE_FINISHED end subroutine get_int subroutine get_real (float) ! get the next real number !--------------------------------------------------------- ! get the next real number !--------------------------------------------------------- Use precision_module implicit none real(dp), intent(out) :: float ! the real number if ( LINE_FINISHED ) then ! assign default of zero float = 0.d0 ; N_CHARS = 1 else ! should be present call locate_next_word ! locate real number read (buffer(START:STOP), *) float ! read float from buffer !K if ( ECHO_KEY ) write (ECHO_UNIT, '(1PE14.6, 1X)', & !K advance='no') float ! append float write (ECHO_UNIT, '(1PE14.6, 1X)', advance='no') float ! append float end if end subroutine get_real END MODULE keyword_buffer ! copyright 1999, J. E. Akin, all rights reserved.