MAKE_DAMPING_MATRIX.f90 File Reference

Go to the source code of this file.

Functions/Subroutines
subroutine	make_damping_matrix (nn, ct, ww, dd, rho, gamma, a11, a12, a13, a21, a22, a23, a31, a32, a33, b11, b12, b13, b21, b22, b23, b31, b32, b33, gg1, gg2, gg3, dd1, dd2, dd3, mc_el, mck_el)
	Make damping matrices.

Function/Subroutine Documentation

make_damping_matrix()

subroutine make_damping_matrix	(	integer*4	nn,
		real*8, dimension(nn)	ct,
		real*8, dimension(nn)	ww,
		real*8, dimension(nn,nn)	dd,
		real*8, dimension(nn,nn,nn)	rho,
		real*8, dimension(nn,nn,nn)	gamma,
		real*8	a11,
		real*8	a12,
		real*8	a13,
		real*8	a21,
		real*8	a22,
		real*8	a23,
		real*8	a31,
		real*8	a32,
		real*8	a33,
		real*8	b11,
		real*8	b12,
		real*8	b13,
		real*8	b21,
		real*8	b22,
		real*8	b23,
		real*8	b31,
		real*8	b32,
		real*8	b33,
		real*8	gg1,
		real*8	gg2,
		real*8	gg3,
		real*8	dd1,
		real*8	dd2,
		real*8	dd3,
		real*8, dimension(nn,nn,nn)	mc_el,
		real*8, dimension(nn,nn,nn)	mck_el
	)

Make damping matrices.

Author

Ilario Mazzieri

Date

September, 2013

Version

1.0

Parameters

[in]	nn	number of 1-D GLL nodes
[in]	xq	GLL nodes
[in]	wq	GLL weights
[in]	dd	spectral derivatives matrix
[in]	rho	mass density
[in]	gamma	damping coefficient
[in]	A11	costant values for the bilinear map
[in]	A12	costant values for the bilinear map
[in]	A13	costant values for the bilinear map
[in]	A21	costant values for the bilinear map
[in]	A22	costant values for the bilinear map
[in]	A23	costant values for the bilinear map
[in]	A31	costant values for the bilinear map
[in]	A32	costant values for the bilinear map
[in]	A33	costant values for the bilinear map
[in]	B11	costant values for the bilinear map
[in]	B12	costant values for the bilinear map
[in]	B13	costant values for the bilinear map
[in]	B21	costant values for the bilinear map
[in]	B22	costant values for the bilinear map
[in]	B23	costant values for the bilinear map
[in]	B31	costant values for the bilinear map
[in]	B32	costant values for the bilinear map
[in]	B33	costant values for the bilinear map
[in]	GG1	costant values for the bilinear map
[in]	GG2	costant values for the bilinear map
[in]	GG3	costant values for the bilinear map
[in]	DD1	costant values for the bilinear map
[in]	DD2	costant values for the bilinear map
[in]	DD3	costant values for the bilinear map
[out]	mc_el	viscous forces proportional to the velocity field
[out]	mck_el	viscous forces proportional to the displacement field

Definition at line 57 of file MAKE_DAMPING_MATRIX.f90.

      
 
      implicit none 
      
      integer*4 :: nn 
      integer*4 :: i,j,k 
 
      real*8 :: a11,a12,a13,a21,a22,a23,a31,a32,a33 
      real*8 :: b11,b12,b13,b21,b22,b23,b31,b32,b33 
      real*8 :: gg1,gg2,gg3,dd1,dd2,dd3 
      real*8 :: dxdx,dxdy,dxdz,dydx,dydy,dydz,dzdx,dzdy,dzdz,det_j 
 
      real*8, dimension(nn) :: ct,ww 
 
      real*8, dimension(nn,nn) :: dd 
 
      real*8, dimension(nn,nn,nn) :: rho,gamma
      real*8, dimension(nn,nn,nn) :: mc_el,mck_el 
 
     
      
      do k = 1,nn 
         do j = 1,nn 
            do i = 1,nn 
               dxdx = a11 +b12*ct(k) +b13*ct(j) & 
                    + gg1*ct(j)*ct(k) 
               dydx = a21 +b22*ct(k) +b23*ct(j) & 
                    + gg2*ct(j)*ct(k) 
               dzdx = a31 +b32*ct(k) +b33*ct(j) & 
                    + gg3*ct(j)*ct(k) 
               
               dxdy = a12 +b11*ct(k) +b13*ct(i) & 
                    + gg1*ct(k)*ct(i) 
               dydy = a22 +b21*ct(k) +b23*ct(i) & 
                    + gg2*ct(k)*ct(i) 
               dzdy = a32 +b31*ct(k) +b33*ct(i) & 
                    + gg3*ct(k)*ct(i) 
               
               dxdz = a13 +b11*ct(j) +b12*ct(i) & 
                    + gg1*ct(i)*ct(j) 
               dydz = a23 +b21*ct(j) +b22*ct(i) & 
                    + gg2*ct(i)*ct(j) 
               dzdz = a33 +b31*ct(j) +b32*ct(i) & 
                    + gg3*ct(i)*ct(j) 
               
               det_j = dxdz * (dydx*dzdy - dzdx*dydy) & 
                     - dydz * (dxdx*dzdy - dzdx*dxdy) & 
                     + dzdz * (dxdx*dydy - dydx*dxdy) 
 
 
             
                    mc_el(i,j,k) = 2 * gamma(i,j,k) * rho(i,j,k) &                                                
                                 * det_j * ww(i) * ww(j) * ww(k) 
                                   
                   mck_el(i,j,k) = (gamma(i,j,k)**2) * rho(i,j,k) & 
                                * det_j * ww(i) * ww(j) * ww(k) 
 
 
            enddo 
         enddo 
      enddo 
      
      return