📄 elmt300.f90
字号:
SUBROUTINE Elmt300( ielem, imats, iswth, iuses )
!.......
! P R O G R A M
!.......
USE CtrlData
USE MeshData
USE ElmtData
USE GlobData
USE ExtraMesh
USE FactorData
IMPLICIT DOUBLE PRECISION ( a-h , o-z )
DIMENSION corep( 3, 8 ), pxgas( 6 ), wxgas( 6 )
DIMENSION smatx( 9, 63 ), pygas( 6 ), wygas( 6 )
DIMENSION dmatx( 6, 6 ), pzgas( 6 ), wzgas( 6 )
DIMENSION snode( 6, 20 ), strgp( 6 )
DIMENSION coreb( 3, 20 ), corgp( 3 )
DIMENSION shape( 4, 21 ), fbody( 3 )
DIMENSION bmatx( 6, 63 ), strng( 6 )
DIMENSION gmatx( 9, 63 ), lnodb( 20 )
DIMENSION ssmat( 9, 9 )
IF( iswth .EQ. 0 ) THEN
CALL RdMat300( imats )
RETURN
ELSE IF( iswth .EQ. 8 ) THEN
nnode = 20
ndofn = 3
RETURN
ELSE IF( iswth .EQ. 9 ) THEN
IF( iuses .EQ. 0 ) RETURN
CALL PlotMesh300( ielem, imats )
RETURN
ELSE IF( iswth .EQ. 21 ) THEN
CALL ExtraMesh300( imats )
RETURN
ELSE IF( iswth .EQ. 22 ) THEN
CALL ExtraDisp300( ielem, imats )
RETURN
END IF
DO idofe = 1, ndofe
DO jdofe = 1, ndofe
stife( idofe, jdofe ) = 0.0D0
END DO
force( idofe ) = 0.0D0
END DO
IF( iuses .EQ. 0 ) RETURN
iflag = 2
nbloc = props( 2, imats ) + 0.5D0
nlinc = props( 3, imats ) + 0.5D0
IF( iswth .EQ. 7 ) WRITE( 12, 2200 ) ielem
IF( iswth .EQ. 1 .OR. iswth .EQ. 5 .OR. iswth .EQ. 10 ) THEN
kelem = lprtx( ielem )
DO ibloc = 1, nbloc
kgaus = 0
igeob = props( ibloc * 2 + 2, imats ) + 0.5D0
imatb = props( ibloc * 2 + 3, imats ) + 0.5D0
ietyp = props( 2, igeob ) + 0.5D0
nxgas = props( 5, igeob ) + 0.5D0
nygas = props( 6, igeob ) + 0.5D0
nzgas = props( 7, igeob ) + 0.5D0
imtyp = props( 16, imatb ) + 0.5D0
CALL Gauss( pxgas, wxgas, nxgas, nerrc )
CALL Gauss( pygas, wygas, nygas, nerrc )
CALL Gauss( pzgas, wzgas, nzgas, nerrc )
IF( nerrc .GT. 0 ) RETURN
CALL Bloc300( coreb, lnodb, igeob, nnodb )
CALL MidPln300( coreb, corep, lnodb, ielem, nnodb )
DO ixgas = 1, nxgas
DO iygas = 1, nygas
DO izgas = 1, nzgas
kgaus = kgaus + 1
xloca = pxgas( ixgas )
yloca = pygas( iygas )
zloca = pzgas( izgas )
xlocd = pxgas( ixgas )
ylocd = pygas( iygas )
CALL Shape300( shape, coreb, lnodb, xloca, yloca, &
zloca, xjaco, nnodb, ielem, iflag )
IF( nerrc .GT. 0 ) RETURN
CALL Bmat300( shape, bmatx, nlinc, 2 )
CALL Dmat300( dmatx, imats, corep, xlocd, ylocd, &
ielem, igeob, imatb, 1 )
IF( nerrc .GT. 0 ) RETURN
dvolu = wxgas( ixgas ) * wygas( iygas ) * &
wzgas( izgas ) * xjaco
IF( iswth .EQ. 1 ) THEN
DO istre = 1, 6
DO idofe = 1, ndofe
smatx( istre, idofe ) = 0.0D0
DO jstre = 1, 6
smatx(istre, idofe) = smatx(istre, idofe) + &
dmatx(istre, jstre) * bmatx(jstre, idofe)
END DO
smatx( istre, idofe ) = smatx( istre, idofe ) * &
dvolu
END DO
END DO
DO idofe = 1, ndofe
DO jdofe = idofe, ndofe
DO istre = 1, 6
stife( idofe, jdofe ) = stife( idofe, jdofe )+ &
bmatx( istre, idofe ) * smatx( istre, jdofe )
END DO
END DO
END DO
IF( nlinc .NE. 0 ) THEN
CALL Gmat300( gmatx, shape, nnode )
CALL Strn300( corgp, strng, shape, nlinc, 3 )
CALL Strs300( corep, strng, strgp, dmatx, xlocd, &
ylocd, imtyp, kgaus, kelem, imatb, &
ietyp, 0 )
CALL SMat300( ssmat, strgp )
IF( nerrc .GT. 0 ) RETURN
DO istre = 1, 9
DO idofe = 1, ndofe
smatx( istre, idofe ) = 0.0D0
DO jstre = 1, 9
smatx(istre, idofe) = smatx(istre, idofe) + &
ssmat(istre, jstre) * gmatx(jstre, idofe)
END DO
smatx(istre, idofe) = smatx(istre, idofe) * &
dvolu
END DO
END DO
DO idofe = 1, ndofe
DO jdofe = idofe, ndofe
DO istre = 1, 9
stife(idofe, jdofe) = stife(idofe, jdofe) + &
gmatx(istre, idofe) * smatx(istre, jdofe)
END DO
END DO
END DO
END IF
END IF
IF( iswth .EQ. 5 .OR. iswth .EQ. 10 .OR. &
nincc .EQ. 2 ) THEN
CALL Strn300( corgp, strng, shape, nlinc, 3 )
CALL Strs300( corep, strng, strgp, dmatx, xlocd, &
ylocd, imtyp, kgaus, kelem, imatb, &
ietyp, 1 )
DO idofe = 1, ndofe
DO istre = 1, 6
force( idofe ) = force( idofe ) - &
bmatx( istre, idofe ) * strgp( istre ) * dvolu
END DO
END DO
END IF
END DO
END DO
END DO
kelem = kelem + 1
END DO
DO idofe = 1, ndofe
DO jdofe = 1, idofe - 1
stife( idofe, jdofe ) = stife( jdofe, idofe )
END DO
END DO
ELSE IF( iswth .EQ. 2 ) THEN
DO ibloc = 1, nbloc
igeob = props( ibloc * 2 + 2, imats ) + 0.5D0
imatb = props( ibloc * 2 + 3, imats ) + 0.5D0
nxgas = props( 5, igeob ) + 0.5D0
nygas = props( 6, igeob ) + 0.5D0
nzgas = props( 7, igeob ) + 0.5D0
imtyp = props( 16, imatb ) + 0.5D0
denst = props( 11, imatb )
CALL Gauss( pxgas, wxgas, nxgas, nerrc )
CALL Gauss( pygas, wygas, nygas, nerrc )
CALL Gauss( pzgas, wzgas, nzgas, nerrc )
IF( nerrc .GT. 0 ) RETURN
CALL Bloc300( coreb, lnodb, igeob, nnodb )
DO ixgas = 1, nxgas
DO iygas = 1, nygas
DO izgas = 1, nzgas
xloca = pxgas( ixgas )
yloca = pygas( iygas )
zloca = pzgas( izgas )
CALL Shape300( shape, coreb, lnodb, xloca, yloca, &
zloca, xjaco, nnodb, ielem, iflag )
IF( nerrc .GT. 0 ) RETURN
dvolu = wxgas( ixgas ) * wygas( iygas ) * &
wzgas( izgas ) * xjaco * denst
DO inode = 1, nnode
iloca = ( inode - 1 ) * 3
DO jnode = 1, nnode
jloca = ( jnode - 1 ) * 3
fwork = shape( 4, inode ) * &
dvolu * shape( 4, jnode )
stife( iloca + 1, jloca + 1 ) = &
stife( iloca + 1, jloca + 1 ) + fwork
stife( iloca + 2, jloca + 2 ) = &
stife( iloca + 2, jloca + 2 ) + fwork
stife( iloca + 3, jloca + 3 ) = &
stife( iloca + 3, jloca + 3 ) + fwork
END DO
END DO
END DO
END DO
END DO
END DO
ELSE IF( iswth .EQ. 3 ) THEN
RETURN
ELSE IF( iswth .EQ. 4 ) THEN
DO ibloc = 1, nbloc
kgaus = 0
igeob = props( ibloc * 2 + 2, imats ) + 0.5D0
imatb = props( ibloc * 2 + 3, imats ) + 0.5D0
ietyp = props( 2, igeob ) + 0.5D0
nxgas = props( 5, igeob ) + 0.5D0
nygas = props( 6, igeob ) + 0.5D0
nzgas = props( 7, igeob ) + 0.5D0
imtyp = props( 16, imatb ) + 0.5D0
CALL Bloc300( coreb, lnodb, igeob, nnodb )
CALL MidPln300( coreb, corep, lnodb, ielem, nnodb )
CALL Gauss( pxgas, wxgas, nxgas, nerrc )
CALL Gauss( pygas, wygas, nygas, nerrc )
CALL Gauss( pzgas, wzgas, nzgas, nerrc )
IF( nerrc .NE. 0 ) RETURN
DO ixgas = 1, nxgas
DO iygas = 1, nygas
DO izgas = 1, nzgas
kgaus = kgaus + 1
xloca = pxgas( ixgas )
yloca = pygas( iygas )
zloca = pzgas( izgas )
xlocd = pxgas( ixgas )
ylocd = pygas( iygas )
CALL Shape300( shape, coreb, lnodb, xloca, yloca, &
zloca, xjaco, nnodb, ielem, iflag )
CALL Dmat300( dmatx, imats, corep, xlocd, ylocd, &
ielem, igeob, imatb, 0 )
IF( nerrc .GT. 0 ) RETURN
dvolu = wxgas( ixgas ) * wygas( iygas ) * &
wzgas( izgas ) * xjaco
CALL Gmat300( gmatx, shape, nnode )
CALL Strn300( corgp, strng, shape, nlinc, 3 )
CALL Strs300( corep, strng, strgp, dmatx, xlocd, &
ylocd, imtyp, kgaus, kelem, imatb, &
ietyp, 0 )
CALL SMat300( ssmat, strgp )
IF( nerrc .GT. 0 ) RETURN
DO istre = 1, 9
DO idofe = 1, ndofe
smatx( istre, idofe ) = 0.0D0
DO jstre = 1, 9
smatx( istre, idofe ) = smatx( istre, idofe ) + &
ssmat( istre, jstre ) * gmatx( jstre, idofe )
END DO
smatx( istre, idofe ) = smatx( istre, idofe ) * &
dvolu
END DO
END DO
DO idofe = 1, ndofe
DO jdofe = idofe, ndofe
DO istre = 1, 9
stife( idofe, jdofe ) = stife( idofe, jdofe ) - &
gmatx( istre, idofe ) * smatx( istre, jdofe )
END DO
END DO
END DO
END DO
END DO
END DO
END DO
DO idofe = 1, ndofe
DO jdofe = 1, idofe - 1
stife( idofe, jdofe ) = stife( jdofe, idofe )
END DO
END DO
RETURN
ELSE IF( MOD( iswth, mswth ) .EQ. 7 ) THEN
kelem = lprtx( ielem )
DO ibloc = 1, nbloc
kgaus = 0
igeob = props( ibloc * 2 + 2, imats ) + 0.5D0
imatb = props( ibloc * 2 + 3, imats ) + 0.5D0
ietyp = props( 2, igeob ) + 0.5D0
nxgas = props( 5, igeob ) + 0.5D0
nygas = props( 6, igeob ) + 0.5D0
nzgas = props( 7, igeob ) + 0.5D0
imtyp = props( 16, imatb ) + 0.5D0
CALL Bloc300( coreb, lnodb, igeob, nnodb )
CALL MidPln300( coreb, corep, lnodb, ielem, nnodb )
CALL Gauss( pxgas, wxgas, nxgas, nerrc )
CALL Gauss( pygas, wygas, nygas, nerrc )
CALL Gauss( pzgas, wzgas, nzgas, nerrc )
volum = 0.0D0
IF( iswth .EQ. mswth + 7 ) THEN
DO inode = 1, nnode
DO jnode = 1, nnode
stife( inode, jnode ) = 0.0D0
END DO
DO istre = 1, 6
snode( istre, inode ) = 0.0D0
smatx( istre, inode ) = 0.0D0
END DO
END DO
END IF
DO ixgas = 1, nxgas
DO iygas = 1, nygas
DO izgas = 1, nzgas
kgaus = kgaus + 1
xloca = pxgas( ixgas )
yloca = pygas( iygas )
zloca = pzgas( izgas )
xlocd = pxgas( ixgas )
ylocd = pygas( iygas )
CALL Shape300( shape, coreb, lnodb, xloca, yloca, &
zloca, xjaco, nnodb, ielem, iflag )
volum = volum + wxgas( ixgas ) * wygas( iygas ) * &
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -