代码
比这篇新的文章:
Demo
比这篇旧的文章: 从NCBI下载基因序列
作者: rui, 点击652次, 评论(0), 收藏者(0), , 打分:
所有评论，共0条：( 我也来说两句)
比这篇旧的文章: 从NCBI下载基因序列
计算Euler方程EigenVectors
语言: Fortran, 标签: Euler 2008/10/30发布 1年前更新作者: rui, 点击652次, 评论(0), 收藏者(0), , 打分:
Fortran语言: 计算Euler方程EigenVectors
001 ! calculate the eigenvectors for \partial f/\partial u
002 subroutine euler_eigen_fu(vx,vy,hm,el,er)
003 implicit none
004 real vx,vy,hm,el(4,4),er(4,4)
005 ! local
006 real q2_a,u_a,v_a,he_a,ssc,gm1
007 real b1,b2
008
009 u_a=vx
010 v_a=vy
011 he_a=hm
012 gm1=eul_gm-1.0
013 q2_a = u_a*u_a + v_a*v_a
014 !c sound speed:
015 ssc = sqrt(gm1*abs(he_a-0.5*q2_a))
016
017 b1 = gm1/(ssc**2)
018 b2 = 0.5*q2_a*b1
019
020 !c right eigenvectors:
021 er(1,1) = 1.0
022 er(1,2) = u_a - ssc
023 er(1,3) = v_a
024 er(1,4) = he_a - u_a*ssc
025 er(2,1) = 0.0
026 er(2,2) = 0.0
027 er(2,3) = 1.0
028 er(2,4) = v_a
029 er(3,1) = 1.0
030 er(3,2) = u_a
031 er(3,3) = v_a
032 er(3,4) = 0.5 * q2_a
033 er(4,1) = 1.0
034 er(4,2) = u_a + ssc
035 er(4,3) = v_a
036 er(4,4) = he_a + u_a*ssc
037 !c left eigenvectors:
038 el(1,1) = 0.5 * (b2+u_a/ssc)
039 el(1,2) = -0.5 * (b1*u_a+1./ssc)
040 el(1,3) = -0.5 * b1 * v_a
041 el(1,4) = 0.5 * b1
042 el(2,1) = -v_a
043 el(2,2) = 0.0
044 el(2,3) = 1.0
045 el(2,4) = 0.0
046 el(3,1) = 1.0 - b2
047 el(3,2) = b1 * u_a
048 el(3,3) = b1 * v_a
049 el(3,4) = -b1
050 el(4,1) = 0.5 * (b2-u_a/ssc)
051 el(4,2) = 0.5 * (-b1*u_a+1.0/ssc)
052 el(4,3) = -0.5 * b1 * v_a
053 el(4,4) = 0.5 * b1
054 !.....................End of Roe Matrix
055 end subroutine euler_eigen_fu
056 ! calculate the eigenvectors for \partial g(u)/\partial u
057 subroutine euler_eigen_gu(vx,vy,hm,el,er)
058 implicit none
059 real vx,vy,hm,el(4,4),er(4,4)
060 ! local
061 real u_a,v_a,he_a,q2_a,gm1
062 real b1,b2,ssc
063
064 u_a=vx
065 v_a=vy
066 he_a=hm
067 gm1=eul_gm-1.0
068
069 q2_a= u_a*u_a + v_a*v_a
070 ! sound speed:
071 ssc = sqrt(gm1*abs(he_a-0.5*q2_a))
072
073 b1 = gm1/(ssc**2)
074 b2 = 0.5*q2_a*b1
075
076 ! right eigenvectors:
077 er(1,1) = 1.0
078 er(1,2) = u_a
079 er(1,3) = v_a - ssc
080 er(1,4) = he_a - v_a*ssc
081 er(2,1) = 0.0
082 er(2,2) = 1.0
083 er(2,3) = 0.0
084 er(2,4) = u_a
085 er(3,1) = 1.0
086 er(3,2) = u_a
087 er(3,3) = v_a
088 er(3,4) = 0.5 * q2_a
089 er(4,1) = 1.0
090 er(4,2) = u_a
091 er(4,3) = v_a + ssc
092 er(4,4) = he_a + v_a*ssc
093 !c left eigenvector:
094 el(1,1) = 0.5 * (b2+v_a/ssc)
095 el(1,2) = -0.5 * b1*u_a
096 el(1,3) = -0.5 * (b1*v_a+1./ssc)
097 el(1,4) = 0.5 * b1
098 el(2,1) = -u_a
099 el(2,2) = 1.0
100 el(2,3) = 0.0
101 el(2,4) = 0.0
102 el(3,1) = 1.0 - b2
103 el(3,2) = b1 * u_a
104 el(3,3) = b1 * v_a
105 el(3,4) = -b1
106 el(4,1) = 0.5 * (b2-v_a/ssc)
107 el(4,2) = -0.5 * b1 * u_a
108 el(4,3) = -0.5 * (b1*v_a-1./ssc)
109 el(4,4) = 0.5 * b1
110 !c.....................End of Roe Matrix...............................
111 end subroutine euler_eigen_gu
002 subroutine euler_eigen_fu(vx,vy,hm,el,er)
003 implicit none
004 real vx,vy,hm,el(4,4),er(4,4)
005 ! local
006 real q2_a,u_a,v_a,he_a,ssc,gm1
007 real b1,b2
008
009 u_a=vx
010 v_a=vy
011 he_a=hm
012 gm1=eul_gm-1.0
013 q2_a = u_a*u_a + v_a*v_a
014 !c sound speed:
015 ssc = sqrt(gm1*abs(he_a-0.5*q2_a))
016
017 b1 = gm1/(ssc**2)
018 b2 = 0.5*q2_a*b1
019
020 !c right eigenvectors:
021 er(1,1) = 1.0
022 er(1,2) = u_a - ssc
023 er(1,3) = v_a
024 er(1,4) = he_a - u_a*ssc
025 er(2,1) = 0.0
026 er(2,2) = 0.0
027 er(2,3) = 1.0
028 er(2,4) = v_a
029 er(3,1) = 1.0
030 er(3,2) = u_a
031 er(3,3) = v_a
032 er(3,4) = 0.5 * q2_a
033 er(4,1) = 1.0
034 er(4,2) = u_a + ssc
035 er(4,3) = v_a
036 er(4,4) = he_a + u_a*ssc
037 !c left eigenvectors:
038 el(1,1) = 0.5 * (b2+u_a/ssc)
039 el(1,2) = -0.5 * (b1*u_a+1./ssc)
040 el(1,3) = -0.5 * b1 * v_a
041 el(1,4) = 0.5 * b1
042 el(2,1) = -v_a
043 el(2,2) = 0.0
044 el(2,3) = 1.0
045 el(2,4) = 0.0
046 el(3,1) = 1.0 - b2
047 el(3,2) = b1 * u_a
048 el(3,3) = b1 * v_a
049 el(3,4) = -b1
050 el(4,1) = 0.5 * (b2-u_a/ssc)
051 el(4,2) = 0.5 * (-b1*u_a+1.0/ssc)
052 el(4,3) = -0.5 * b1 * v_a
053 el(4,4) = 0.5 * b1
054 !.....................End of Roe Matrix
055 end subroutine euler_eigen_fu
056 ! calculate the eigenvectors for \partial g(u)/\partial u
057 subroutine euler_eigen_gu(vx,vy,hm,el,er)
058 implicit none
059 real vx,vy,hm,el(4,4),er(4,4)
060 ! local
061 real u_a,v_a,he_a,q2_a,gm1
062 real b1,b2,ssc
063
064 u_a=vx
065 v_a=vy
066 he_a=hm
067 gm1=eul_gm-1.0
068
069 q2_a= u_a*u_a + v_a*v_a
070 ! sound speed:
071 ssc = sqrt(gm1*abs(he_a-0.5*q2_a))
072
073 b1 = gm1/(ssc**2)
074 b2 = 0.5*q2_a*b1
075
076 ! right eigenvectors:
077 er(1,1) = 1.0
078 er(1,2) = u_a
079 er(1,3) = v_a - ssc
080 er(1,4) = he_a - v_a*ssc
081 er(2,1) = 0.0
082 er(2,2) = 1.0
083 er(2,3) = 0.0
084 er(2,4) = u_a
085 er(3,1) = 1.0
086 er(3,2) = u_a
087 er(3,3) = v_a
088 er(3,4) = 0.5 * q2_a
089 er(4,1) = 1.0
090 er(4,2) = u_a
091 er(4,3) = v_a + ssc
092 er(4,4) = he_a + v_a*ssc
093 !c left eigenvector:
094 el(1,1) = 0.5 * (b2+v_a/ssc)
095 el(1,2) = -0.5 * b1*u_a
096 el(1,3) = -0.5 * (b1*v_a+1./ssc)
097 el(1,4) = 0.5 * b1
098 el(2,1) = -u_a
099 el(2,2) = 1.0
100 el(2,3) = 0.0
101 el(2,4) = 0.0
102 el(3,1) = 1.0 - b2
103 el(3,2) = b1 * u_a
104 el(3,3) = b1 * v_a
105 el(3,4) = -b1
106 el(4,1) = 0.5 * (b2-v_a/ssc)
107 el(4,2) = -0.5 * b1 * u_a
108 el(4,3) = -0.5 * (b1*v_a-1./ssc)
109 el(4,4) = 0.5 * b1
110 !c.....................End of Roe Matrix...............................
111 end subroutine euler_eigen_gu
所有评论，共0条：( 我也来说两句)