悲催的科学匠人 - 冷水's blog
Latex Listings Package 使用方法[ZZ]
转自 http://xavier.yo2.cn/articles/latex-listings-package-%E4%BD%BF%E7%94%A8%E6%96%B9%E6%B3%95.html
这里就不班门弄斧了，直接转载过来备忘好了。
转载自http://blog.linuxgem.org/lyanry/show/319.html
listings 是专用于代码排版的 LaTeX 宏包，可对关键词、注释和字符串等使用不同的字体和颜色或颜色，也可以为代码添加边框、背景等风格。
1 基本用法
下面给出一份用于排版 C 语言 HelloWorld 程序代码的完整的 LaTeX 文档：
\documentclass{article}
\usepackage{listings}
\begin{document}
\begin{lstlisting}[language=C]
int main(int argc, char ** argv)
{
printf("Hello world!\n");
return 0;
}
\end{lstlisting}
\end{document}
注意，要使用 listings 宏包提供的语法高亮，需要 xcolor 宏包支持。
4 添加边框
listings 宏包为代码边框提供了很多风格，大体可分为带有阴影的边框与圆角边框。这里仅仅给出一个阴影边框的示例，至于其它边框风格，可查阅 listings 宏包文档，里面给出了一些示例。
下面 LaTeX 源文档将为代码添加阴影边框，并将阴影设置为浅灰色：
\begin{lstlisting}[language={[ANSI]C}, keywordstyle=\color{blue!70}, commentstyle=\color{red!50!green!50!blue!50}, frame=shadowbox, rulesepcolor=\color{red!20!green!20!blue!20}]
int main(int argc, char ** argv)
{
printf("Hello world!\n") ;
return 0;
}
\end{lstlisting}
5 添加行号
很多时候需要对文档中的代码进行解释，只有带有行号的代码才可以让解释更清晰，因为你只需要说第 x 行代码有什么作用即可。如果没有行号，那对读者而言就太残忍了，他们不得不从你的文字叙述中得知行号信息，然后去一行一行的查到相应代码行。
listings 宏包通过参数 numbers 来设定行号，该参数的值有两个，分别是 left 与 right，表示行号显示在代码的左侧还是右侧。下面为带有边框的代码添加行号，并设置行号字体为 \tiny：
\begin{lstlisting}[language={[ANSI]C}, numbers=left, numberstyle=\tiny, keywordstyle=\color{blue!70}, commentstyle=\color{red!50!green!50!blue!50}, frame=shadowbox, rulesepcolor=\color{red!20!green!20!blue!20}]
int main(int argc, char ** argv)
{
printf("Hello world!\n") ;
return 0;
}
\end{lstlisting}
6 全局设置
上面所给的各个示例中，lstlisting 环境后面尾随了很多参数，要是每使用一次 lstlisting 环境就要设置这么多参数，那就没什么意思了。
可以使用 \lstset 命令在 LaTeX 源文档的导言区设定好 lstlisting 环境所用的公共参数，如下：
\documentclass{article}
\usepackage{listings}
\usepackage{xcolor}
\begin{document}
\lstset{numbers=left,
numberstyle=\tiny,
keywordstyle=\color{blue!70}, commentstyle=\color{red!50!green!50!blue!50},
frame=shadowbox,
rulesepcolor=\color{red!20!green!20!blue!20}
}
\begin{lstlisting}[language={[ANSI]C}]
int main(int argc, char ** argv)
{
printf("Hello world!\n") ;
return 0;
}
\end{lstlisting}
\end{document}
7 显示中文
listings 宏包默认是不支持包含中文字串的代码显示的，但是可以使用 “逃逸” 字串来显示中文。
在 \lstset 命令中设置逃逸字串的开始符号与终止符号，推荐使用的符号是左引号，即 “ `”
\lstset{numbers=left,
numberstyle=\tiny, keywordstyle=\color{blue!70}, commentstyle=\color{red!50!green!50!blue!50},
frame=shadowbox, rulesepcolor=\color{red!20!green!20!blue!20},
escapeinside=``}
……
\begin{lstlisting}[language={[ANSI]C}]
int main(int argc, char ** argv)
{
printf("`我爱中文`!\n") ;
return 0;
}
\end{lstlisting}
8 调整一下边距
listings 的代码框的宽度默认是与页芯等宽的，其上边距也过于小，可根据自己的审美观念适度调整一下。我通常是将代码框的左右边距设置为 2em，上边距为 1em，下边距采用默认值即可，所作设定如下：
\lstset{numbers=left, numberstyle=\tiny,
keywordstyle=\color{blue!70},
commentstyle=\color{red!50!green!50!blue!50},
frame=shadowbox,
rulesepcolor=\color{red!20!green!20!blue!20},escapeinside=``,
xleftmargin=2em,xrightmargin=2em, aboveskip=1em}
PLOT3D文件格式[ZZ]
http://hi.baidu.com/necrohan/blog/item/4997711f294102fde0fe0b8e.html
其实早已熟悉，只是借来存档一下。此外，以前不知iblank的unformatted格式写法，这个博文居然有。
PLOT3D文件格式
FLOW3D可以导出的数据格式有限，说明中说支持tecplot，但是是需要tecplot插件的，而且要360版本才行，找了2个插件可是不能用，而且插件的功能好像很有限。
PLOT3D格式的数据可以在后处理时导出，为此查了一下PLOT3D格式的数据，实在不行再转为tecplot格式的数据。
PLOT3D 数据格式
PLOT3D数据格式源于NASA，广泛用于规则网格的CFD数据文件。PLOT3D文件可以是ASCII, 也可是Fortran unformatted 或 C binary形式。
PLOT3D文件分为网格文件(XYZ 文件), 空气动力学结果文件 (Q 文件)和通用结果文件(函数文件 + 函数名称文件)。网格文件中可加入所谓的IBlank参数。
IBlank参数的定义：
IBlank是在每一网格点上给出的一个正数值，定义如下：
0 - 计算域之外，非流体点
1 - 正常点
2 - 固面边界点
负数 - 分块网格界面点，其数值为相邻网格块的块号
以下为各文件使用FORTRAN读入的语句，所有文件均为无格式文件。
网格文件（XYZ文件）：
XYZ 文件, 单块（single-block）:
READ(IUNIT) IMAX,JMAX,KMAX
READ(IUNIT) (((X(I,J,K),I=1,IMAX),J=1,JMAX),K=1,KMAX),&
(((Y(I,J,K),I=1,IMAX),J=1,JMAX),K=1,KMAX),&
(((Z(I,J,K),I=1,IMAX),J=1,JMAX),K=1,KMAX)
XYZ 文件, 单块（single-block）, 加 IBlank:
READ(IUNIT) IMAX,JMAX,KMAX
READ(IUNIT) (((X(I,J,K),I=1,IMAX),J=1,JMAX),K=1,KMAX),&
(((Y(I,J,K),I=1,IMAX),J=1,JMAX),K=1,KMAX),&
(((Z(I,J,K),I=1,IMAX),J=1,JMAX),K=1,KMAX),&
(((IBLANK(I,J,K),I=1,IMAX),J=1,JMAX),K=1,KMAX)
XYZ 文件, 二维, 单块（single-block）:
READ(IUNIT) IMAX,JMAX
READ(IUNIT) ((X(I,J),I=1,IMAX),J=1,JMAX),&
((Y(I,J),I=1,IMAX),J=1,JMAX)
XYZ 文件, 多块（multi-block）
READ(IUNIT) NBLOCK
READ(IUNIT) (IMAX(N),JMAX(N),KMAX(N),N=1,NBLOCK)
DO N=1,NBLOCK
READ(IUNIT) (((X(I,J,K),I=1,IMAX(N)),J=1,JMAX(N)),K=1,KMAX(N)),&
(((Y(I,J,K),I=1,IMAX(N)),J=1,JMAX(N)),K=1,KMAX(N)),&
(((Z(I,J,K),I=1,IMAX(N)),J=1,JMAX(N)),K=1,KMAX(N))
ENDDO
XYZ 文件, 多块（multi-block）, 加 IBlank:
READ(IUNIT) NBLOCK
READ(IUNIT) (IMAX(N),JMAX(N),KMAX(N),N=1,NBLOCK)
DO N=1,NBLOCK
READ(IUNIT) (((X(I,J,K),I=1,IMAX(N)),J=1,JMAX(N)),K=1,KMAX(N)),&
(((Y(I,J,K),I=1,IMAX(N)),J=1,JMAX(N)),K=1,KMAX(N)),&
(((Z(I,J,K),I=1,IMAX(N)),J=1,JMAX(N)),K=1,KMAX(N)),&
(((IBLANK(I,J,K),I=1,IMAX(N)),J=1,JMAX(N)),K=1,KMAX(N))
ENDDO
XYZ 文件, 二维, 多块（multi-block）
READ(IUNIT) NBLOCK
READ(IUNIT) (IMAX(N),JMAX(N),N=1,NBLOCK)
DO N=1,NBLOCK
READ(IUNIT) ((X(I,J),I=1,IMAX(N)),J=1,JMAX(N)),&
((Y(I,J),I=1,IMAX(N)),J=1,JMAX(N))
ENDDO
空气动力学结果文件 (Q 文件)(Q 文件定义过窄，现已很少使用。):
Q 文件专为外流空气动力学设计，对三维流动，数组变量如下：
Q1 - 无量纲 密度
Q2 - 无量纲 X-动量
Q3 - 无量纲 Y-动量
Q4 - 无量纲 Z-动量
Q5 - 无量纲 总能
另加4个参数：
FSMACH - 自由流马赫数
ALPHA - 攻角
RE - 雷诺数
TIME - 时间
Q 文件, 单块（single-block）:
READ(IUNIT) IMAX,JMAX,KMAX
READ(IUNIT) FSMACH,ALPHA,RE,TIME
READ(IUNIT) ((((Q(I,J,K,M),I=1,IMAX),J=1,JMAX),K=1,KMAX),M=1,5)
Q 文件, 二维, 单块（single-block）:
READ(IUNIT) IMAX,JMAX
READ(IUNIT) FSMACH,ALPHA,RE,TIME
READ(IUNIT) (((Q(I,J,M),I=1,IMAX),J=1,JMAX),M=1,4)
Q 文件, 多块（multi-block）:
READ(IUNIT) NBLOCK
READ(IUNIT) (IMAX(N),JMAX(N),KMAX(N),N=1,NBLOCK)
DO N=1,NBLOCK
READ(IUNIT) FSMACH,ALPHA,RE,TIME
READ(IUNIT)((((Q(I,J,K,M),I=1,IMAX(N)),J=1,JMAX(N)),K=1,KMAX(N)),M=1,5)
ENDDO
Q 文件定义过窄，现已很少使用。大多CFD工作者使用通用结果文件，即
所谓的函数文件（function file）。函数文件中可定义任意数量的数组变量，
其名称在另一函数名称文件（function name file）中给出。
函数文件:
函数文件, 单块（single-block）:
READ(IUNIT) IMAX,JMAX,KMAX
READ(IUNIT) ((((F(I,J,K,M),I=1,IMAX),J=1,JMAX),K=1,KMAX),M=1,NFUN)
函数文件, 二维, 单块（single-block）:
READ(IUNIT) IMAX,JMAX
READ(IUNIT) (((F(I,J,M),I=1,IMAX),J=1,JMAX),M=1,NFUN)
函数文件, 多块（multi-block）:
READ(IUNIT) NBLOCK
READ(IUNIT) (IMAX(N),JMAX(N),KMAX(N),N=1,NBLOCK)
DO N=1,NBLOCK
READ(IUNIT) ((((F(I,J,K,M),I=1,IMAX(N)),J=1,JMAX(N)),K=1,KMAX(N)),M=1,NFUN)
ENDDO
函数名称文件:
函数名称文件是一ASCII文件，列有函数文件中数组变量的对应名称。以下为一例子：
density
pressure
u;velocity vector
v
w
temperature
turbulence energy
...
...
注意事项：
1 名称的数量和排列次序与函数文件中数组变量相同。
2 分割号 ";" 表示向量行的开始。";"之右为向量名，";"之左为X分量。其下两行为Y分量和Z分量。
http://www.grc.nasa.gov/WWW/wind/valid/plot3d.html
下面例子假设网格文件句柄为7，结果文件句柄为8。
2D, Whole, Formatted, Single-Block Grid and Solution
parameter ( imax = 100 )
parameter ( jmax = 100 )
integer i
integer j
integer m
integer n
integer ni
integer nj
real mach ! freestream Mach number
real alpha ! freestream angle-of-attack
real reyn ! freestream Reynolds number
real time ! time
real x(imax,jmax)
real y(imax,jmax)
real q(imax,jmax,4)
open ( unit=7, form='formatted', file='2D.x' )
open ( unit=8, form='formatted', file='2D.q' )
read(7,*) ni, nj
read(7,*)
& (( x(i,j), i=1,ni), j=1,nj),
& (( y(i,j), i=1,ni), j=1,nj)
read(8,*) ni, nj
read(8,*) mach, alpha, reyn, time
read(8,*) ((( q(i,j,n), i=1,ni), j=1,nj), n=1,4)
3D, Whole, Unformatted, Multi-Block Grid and Solution
parameter ( imax = 100 )
parameter ( jmax = 100 )
parameter ( kmax = 100 )
parameter ( nbmax = 10 )
integer i
integer j
integer m
integer n
integer nblocks
integer ni (nbmax)
integer nj (nbmax)
integer nk (nbmax)
real mach ! freestream Mach number
real alpha ! freestream angle-of-attack
real reyn ! freestream Reynolds number
real time ! time
real x(imax,jmax,kmax,nbmax)
real y(imax,jmax,kmax,nbmax)
real z(imax,jmax,kmax,nbmax)
real q(imax,jmax,kmax,nbmax,5)
open ( unit=7, form='unformatted', file='3D.x' )
open ( unit=8, form='unformatted', file='3D.q' )
read(7) nblocks
read(7) ( ni(m), nj(m), nk(m), m = 1, nblocks )
do m = 1, nblocks
read(7)
& ((( x(i,j,k,m), i=1,ni(m)), j=1,nj(m)), k=1,nk(m)),
& ((( y(i,j,k,m), i=1,ni(m)), j=1,nj(m)), k=1,nk(m)),
& ((( z(i,j,k,m), i=1,ni(m)), j=1,nj(m)), k=1,nk(m))
enddo
read(8) nblocks
read(8) ( ni(m), nj(m), nk(m), m = 1, nblocks )
do m = 1, nblocks
read(8) mach, alpha, reyn, time
read(8)
& (((( q(i,j,k,m,n), i=1,ni(m)), j=1,nj(m)), k=1,nk(m)), n=1,5)
enddo
TECPLOT 的 PREPLOT 可以把 PLOT3D 文件转为 TECPLOT 格式。（没有测试）
20090526 add
由于后来发现可以用flow3d自带的命令行方式后处理提取数据，所以这个方法放弃。
PLOT3D格式的数据可以在后处理时导出，为此查了一下PLOT3D格式的数据，实在不行再转为tecplot格式的数据。
PLOT3D 数据格式
PLOT3D数据格式源于NASA，广泛用于规则网格的CFD数据文件。PLOT3D文件可以是ASCII, 也可是Fortran unformatted 或 C binary形式。
PLOT3D文件分为网格文件(XYZ 文件), 空气动力学结果文件 (Q 文件)和通用结果文件(函数文件 + 函数名称文件)。网格文件中可加入所谓的IBlank参数。
IBlank参数的定义：
IBlank是在每一网格点上给出的一个正数值，定义如下：
0 - 计算域之外，非流体点
1 - 正常点
2 - 固面边界点
负数 - 分块网格界面点，其数值为相邻网格块的块号
以下为各文件使用FORTRAN读入的语句，所有文件均为无格式文件。
网格文件（XYZ文件）：
XYZ 文件, 单块（single-block）:
READ(IUNIT) IMAX,JMAX,KMAX
READ(IUNIT) (((X(I,J,K),I=1,IMAX),J=1,JMAX),K=1,KMAX),&
(((Y(I,J,K),I=1,IMAX),J=1,JMAX),K=1,KMAX),&
(((Z(I,J,K),I=1,IMAX),J=1,JMAX),K=1,KMAX)
XYZ 文件, 单块（single-block）, 加 IBlank:
READ(IUNIT) IMAX,JMAX,KMAX
READ(IUNIT) (((X(I,J,K),I=1,IMAX),J=1,JMAX),K=1,KMAX),&
(((Y(I,J,K),I=1,IMAX),J=1,JMAX),K=1,KMAX),&
(((Z(I,J,K),I=1,IMAX),J=1,JMAX),K=1,KMAX),&
(((IBLANK(I,J,K),I=1,IMAX),J=1,JMAX),K=1,KMAX)
XYZ 文件, 二维, 单块（single-block）:
READ(IUNIT) IMAX,JMAX
READ(IUNIT) ((X(I,J),I=1,IMAX),J=1,JMAX),&
((Y(I,J),I=1,IMAX),J=1,JMAX)
XYZ 文件, 多块（multi-block）
READ(IUNIT) NBLOCK
READ(IUNIT) (IMAX(N),JMAX(N),KMAX(N),N=1,NBLOCK)
DO N=1,NBLOCK
READ(IUNIT) (((X(I,J,K),I=1,IMAX(N)),J=1,JMAX(N)),K=1,KMAX(N)),&
(((Y(I,J,K),I=1,IMAX(N)),J=1,JMAX(N)),K=1,KMAX(N)),&
(((Z(I,J,K),I=1,IMAX(N)),J=1,JMAX(N)),K=1,KMAX(N))
ENDDO
XYZ 文件, 多块（multi-block）, 加 IBlank:
READ(IUNIT) NBLOCK
READ(IUNIT) (IMAX(N),JMAX(N),KMAX(N),N=1,NBLOCK)
DO N=1,NBLOCK
READ(IUNIT) (((X(I,J,K),I=1,IMAX(N)),J=1,JMAX(N)),K=1,KMAX(N)),&
(((Y(I,J,K),I=1,IMAX(N)),J=1,JMAX(N)),K=1,KMAX(N)),&
(((Z(I,J,K),I=1,IMAX(N)),J=1,JMAX(N)),K=1,KMAX(N)),&
(((IBLANK(I,J,K),I=1,IMAX(N)),J=1,JMAX(N)),K=1,KMAX(N))
ENDDO
XYZ 文件, 二维, 多块（multi-block）
READ(IUNIT) NBLOCK
READ(IUNIT) (IMAX(N),JMAX(N),N=1,NBLOCK)
DO N=1,NBLOCK
READ(IUNIT) ((X(I,J),I=1,IMAX(N)),J=1,JMAX(N)),&
((Y(I,J),I=1,IMAX(N)),J=1,JMAX(N))
ENDDO
空气动力学结果文件 (Q 文件)(Q 文件定义过窄，现已很少使用。):
Q 文件专为外流空气动力学设计，对三维流动，数组变量如下：
Q1 - 无量纲 密度
Q2 - 无量纲 X-动量
Q3 - 无量纲 Y-动量
Q4 - 无量纲 Z-动量
Q5 - 无量纲 总能
另加4个参数：
FSMACH - 自由流马赫数
ALPHA - 攻角
RE - 雷诺数
TIME - 时间
Q 文件, 单块（single-block）:
READ(IUNIT) IMAX,JMAX,KMAX
READ(IUNIT) FSMACH,ALPHA,RE,TIME
READ(IUNIT) ((((Q(I,J,K,M),I=1,IMAX),J=1,JMAX),K=1,KMAX),M=1,5)
Q 文件, 二维, 单块（single-block）:
READ(IUNIT) IMAX,JMAX
READ(IUNIT) FSMACH,ALPHA,RE,TIME
READ(IUNIT) (((Q(I,J,M),I=1,IMAX),J=1,JMAX),M=1,4)
Q 文件, 多块（multi-block）:
READ(IUNIT) NBLOCK
READ(IUNIT) (IMAX(N),JMAX(N),KMAX(N),N=1,NBLOCK)
DO N=1,NBLOCK
READ(IUNIT) FSMACH,ALPHA,RE,TIME
READ(IUNIT)((((Q(I,J,K,M),I=1,IMAX(N)),J=1,JMAX(N)),K=1,KMAX(N)),M=1,5)
ENDDO
Q 文件定义过窄，现已很少使用。大多CFD工作者使用通用结果文件，即
所谓的函数文件（function file）。函数文件中可定义任意数量的数组变量，
其名称在另一函数名称文件（function name file）中给出。
函数文件:
函数文件, 单块（single-block）:
READ(IUNIT) IMAX,JMAX,KMAX
READ(IUNIT) ((((F(I,J,K,M),I=1,IMAX),J=1,JMAX),K=1,KMAX),M=1,NFUN)
函数文件, 二维, 单块（single-block）:
READ(IUNIT) IMAX,JMAX
READ(IUNIT) (((F(I,J,M),I=1,IMAX),J=1,JMAX),M=1,NFUN)
函数文件, 多块（multi-block）:
READ(IUNIT) NBLOCK
READ(IUNIT) (IMAX(N),JMAX(N),KMAX(N),N=1,NBLOCK)
DO N=1,NBLOCK
READ(IUNIT) ((((F(I,J,K,M),I=1,IMAX(N)),J=1,JMAX(N)),K=1,KMAX(N)),M=1,NFUN)
ENDDO
函数名称文件:
函数名称文件是一ASCII文件，列有函数文件中数组变量的对应名称。以下为一例子：
density
pressure
u;velocity vector
v
w
temperature
turbulence energy
...
...
注意事项：
1 名称的数量和排列次序与函数文件中数组变量相同。
2 分割号 ";" 表示向量行的开始。";"之右为向量名，";"之左为X分量。其下两行为Y分量和Z分量。
http://www.grc.nasa.gov/WWW/wind/valid/plot3d.html
下面例子假设网格文件句柄为7，结果文件句柄为8。
2D, Whole, Formatted, Single-Block Grid and Solution
parameter ( imax = 100 )
parameter ( jmax = 100 )
integer i
integer j
integer m
integer n
integer ni
integer nj
real mach ! freestream Mach number
real alpha ! freestream angle-of-attack
real reyn ! freestream Reynolds number
real time ! time
real x(imax,jmax)
real y(imax,jmax)
real q(imax,jmax,4)
open ( unit=7, form='formatted', file='2D.x' )
open ( unit=8, form='formatted', file='2D.q' )
read(7,*) ni, nj
read(7,*)
& (( x(i,j), i=1,ni), j=1,nj),
& (( y(i,j), i=1,ni), j=1,nj)
read(8,*) ni, nj
read(8,*) mach, alpha, reyn, time
read(8,*) ((( q(i,j,n), i=1,ni), j=1,nj), n=1,4)
3D, Whole, Unformatted, Multi-Block Grid and Solution
parameter ( imax = 100 )
parameter ( jmax = 100 )
parameter ( kmax = 100 )
parameter ( nbmax = 10 )
integer i
integer j
integer m
integer n
integer nblocks
integer ni (nbmax)
integer nj (nbmax)
integer nk (nbmax)
real mach ! freestream Mach number
real alpha ! freestream angle-of-attack
real reyn ! freestream Reynolds number
real time ! time
real x(imax,jmax,kmax,nbmax)
real y(imax,jmax,kmax,nbmax)
real z(imax,jmax,kmax,nbmax)
real q(imax,jmax,kmax,nbmax,5)
open ( unit=7, form='unformatted', file='3D.x' )
open ( unit=8, form='unformatted', file='3D.q' )
read(7) nblocks
read(7) ( ni(m), nj(m), nk(m), m = 1, nblocks )
do m = 1, nblocks
read(7)
& ((( x(i,j,k,m), i=1,ni(m)), j=1,nj(m)), k=1,nk(m)),
& ((( y(i,j,k,m), i=1,ni(m)), j=1,nj(m)), k=1,nk(m)),
& ((( z(i,j,k,m), i=1,ni(m)), j=1,nj(m)), k=1,nk(m))
enddo
read(8) nblocks
read(8) ( ni(m), nj(m), nk(m), m = 1, nblocks )
do m = 1, nblocks
read(8) mach, alpha, reyn, time
read(8)
& (((( q(i,j,k,m,n), i=1,ni(m)), j=1,nj(m)), k=1,nk(m)), n=1,5)
enddo
TECPLOT 的 PREPLOT 可以把 PLOT3D 文件转为 TECPLOT 格式。（没有测试）
20090526 add
由于后来发现可以用flow3d自带的命令行方式后处理提取数据，所以这个方法放弃。
一个简单的动态数据容器
能够动态的添加成员数据，目前只支持integer(4)和real(8)
! dynamic data type
! Qu Kun, April 2012
! 一个简单的动态数据结构容器，能够容纳INTEGER(4)和REAL(8)数据，维数范围为0到5
! 多个容器可链接成为链表
! 在链表中创建插入容器 dyd_new，指定数据类型、维数、指标范围和名称
! 获取容器指针 dyd_findnode 在链表中搜索指定名称的容器
! 释放容器内的数据空间 dyd_empty
! 开辟容器数据空间 dyd_construct
! 获取数据指针 dyd_getdata， 通过名称得到数据指针
!
! 一般使用最常用的是 dyd_new和dyd_getdata
! dyd_new的使用方法，建议采用关键字方式指定参数
! 1 变量指定数组范围，range11,range12,range21,range22,...,range51,range52
! CALL dyd_new(head=list, elementtype='integer(4)', dataname='X', ierr=ierr, &
! range11=0,range12=9)
! 2 字符串指定数组范围
! CALL dyd_new(head=list, elementtype='integer(4)', arrayrange='(0:9)', &
! dataname='X', ierr=ierr)
! 3 数组方式指定数组范围
! CALL dyd_new(head=list, atype=int4, adim=1,arange=irange, aname='X')
! 比较方便的是前2种
!
! 获取数据指针的方法为
! INTEGER(4),DIMENSION(:),POINTER :: X
! ...
! CALL getdata(list,'X', X)
! DO i=0,9; X(i) = i; ENDDO
MODULE qkdynadata
IMPLICIT NONE
INTEGER,PARAMETER,PRIVATE :: dyd_maxdim=5,dyd_namelen=32
INTEGER,PARAMETER,PRIVATE :: int4=1,real8=2
TYPE dynadata
INTEGER :: mdim=-1
INTEGER :: mrange(2,dyd_maxdim)=-1
INTEGER :: metype=-1 ! data type of the element, integer(4), real(8),real(4)
CHARACTER(LEN=dyd_namelen) mname
INTEGER(4),POINTER :: mint4d0p=>NULL()
INTEGER(4),POINTER,DIMENSION(:) :: mint4d1p=>NULL()
INTEGER(4),POINTER,DIMENSION(:,:) :: mint4d2p=>NULL()
INTEGER(4),POINTER,DIMENSION(:,:,:) :: mint4d3p=>NULL()
INTEGER(4),POINTER,DIMENSION(:,:,:,:) :: mint4d4p=>NULL()
INTEGER(4),POINTER,DIMENSION(:,:,:,:,:) :: mint4d5p=>NULL()
REAL(8),POINTER :: mreal8d0p=>NULL()
REAL(8),POINTER,DIMENSION(:) :: mreal8d1p=>NULL()
REAL(8),POINTER,DIMENSION(:,:) :: mreal8d2p=>NULL()
REAL(8),POINTER,DIMENSION(:,:,:) :: mreal8d3p=>NULL()
REAL(8),POINTER,DIMENSION(:,:,:,:) :: mreal8d4p=>NULL()
REAL(8),POINTER,DIMENSION(:,:,:,:,:) :: mreal8d5p=>NULL()
!pointer for inner dynamic data member
TYPE(dynadata),POINTER :: mdynadata
! pointer for linklist
TYPE(dynadata),POINTER :: mprev,mnext,mhead
END TYPE
INTERFACE getdata
MODULE PROCEDURE getdata_int4p0
MODULE PROCEDURE getdata_int4p1
MODULE PROCEDURE getdata_int4p2
MODULE PROCEDURE getdata_int4p3
MODULE PROCEDURE getdata_int4p4
MODULE PROCEDURE getdata_int4p5
MODULE PROCEDURE getdata_real8p0
MODULE PROCEDURE getdata_real8p1
MODULE PROCEDURE getdata_real8p2
MODULE PROCEDURE getdata_real8p3
MODULE PROCEDURE getdata_real8p4
MODULE PROCEDURE getdata_real8p5
END INTERFACE
CONTAINS
SUBROUTINE getdata_int4p0(head,aname, dp)
INTEGER(4),POINTER :: dp
TYPE(dynadata),POINTER,INTENT(IN) :: head
CHARACTER(LEN=*),INTENT(IN) :: aname
!
TYPE(dynadata),POINTER :: iter=>NULL()
NULLIFY(dp)
iter => dyd_findnode(head,aname)
IF(ASSOCIATED(iter)) dp => iter%mint4d0p
END SUBROUTINE
SUBROUTINE getdata_int4p1(head,aname,dp)
INTEGER(4),DIMENSION(:),POINTER ::dp
TYPE(dynadata),POINTER,INTENT(IN) :: head
CHARACTER(LEN=*),INTENT(IN) :: aname
!
TYPE(dynadata),POINTER :: iter=>NULL()
NULLIFY(dp)
iter => dyd_findnode(head,aname)
! WRITE(*,*) iter%mname
dp => iter%mint4d1p
END SUBROUTINE
SUBROUTINE getdata_int4p2(head,aname,dp)
INTEGER(4),DIMENSION(:,:),POINTER ::dp
TYPE(dynadata),POINTER,INTENT(IN) :: head
CHARACTER(LEN=*),INTENT(IN) :: aname
!
TYPE(dynadata),POINTER :: iter=>NULL()
NULLIFY(dp)
iter => dyd_findnode(head,aname)
dp => iter%mint4d2p
END SUBROUTINE
SUBROUTINE getdata_int4p3(head,aname,dp)
INTEGER(4),DIMENSION(:,:,:),POINTER ::dp
TYPE(dynadata),POINTER,INTENT(IN) :: head
CHARACTER(LEN=*),INTENT(IN) :: aname
!
TYPE(dynadata),POINTER :: iter=>NULL()
NULLIFY(dp)
iter => dyd_findnode(head,aname)
dp => iter%mint4d3p
END SUBROUTINE
SUBROUTINE getdata_int4p4(head,aname,dp)
INTEGER(4),DIMENSION(:,:,:,:),POINTER ::dp
TYPE(dynadata),POINTER,INTENT(IN) :: head
CHARACTER(LEN=*),INTENT(IN) :: aname
!
TYPE(dynadata),POINTER :: iter=>NULL()
NULLIFY(dp)
iter => dyd_findnode(head,aname)
dp => iter%mint4d4p
END SUBROUTINE
SUBROUTINE getdata_int4p5(head,aname,dp)
INTEGER(4),DIMENSION(:,:,:,:,:),POINTER ::dp
TYPE(dynadata),POINTER,INTENT(IN) :: head
CHARACTER(LEN=*),INTENT(IN) :: aname
!
TYPE(dynadata),POINTER :: iter=>NULL()
NULLIFY(dp)
iter => dyd_findnode(head,aname)
dp => iter%mint4d5p
END SUBROUTINE
SUBROUTINE getdata_real8p0(head,aname, dp)
REAL(8),POINTER :: dp
TYPE(dynadata),POINTER,INTENT(IN) :: head
CHARACTER(LEN=*),INTENT(IN) :: aname
!
TYPE(dynadata),POINTER :: iter=>NULL()
NULLIFY(dp)
iter => dyd_findnode(head,aname)
dp => iter%mreal8d0p
END SUBROUTINE
SUBROUTINE getdata_real8p1(head,aname,dp)
REAL(8),DIMENSION(:),POINTER ::dp
TYPE(dynadata),POINTER,INTENT(IN) :: head
CHARACTER(LEN=*),INTENT(IN) :: aname
!
TYPE(dynadata),POINTER :: iter=>NULL()
NULLIFY(dp)
iter => dyd_findnode(head,aname)
dp => iter%mreal8d1p
END SUBROUTINE
SUBROUTINE getdata_real8p2(head,aname,dp)
REAL(8),DIMENSION(:,:),POINTER ::dp
TYPE(dynadata),POINTER,INTENT(IN) :: head
CHARACTER(LEN=*),INTENT(IN) :: aname
!
TYPE(dynadata),POINTER :: iter=>NULL()
NULLIFY(dp)
iter => dyd_findnode(head,aname)
dp => iter%mreal8d2p
END SUBROUTINE
SUBROUTINE getdata_real8p3(head,aname,dp)
REAL(8),DIMENSION(:,:,:),POINTER ::dp
TYPE(dynadata),POINTER,INTENT(IN) :: head
CHARACTER(LEN=*),INTENT(IN) :: aname
!
TYPE(dynadata),POINTER :: iter=>NULL()
NULLIFY(dp)
iter => dyd_findnode(head,aname)
dp => iter%mreal8d3p
END SUBROUTINE
SUBROUTINE getdata_real8p4(head,aname,dp)
REAL(8),DIMENSION(:,:,:,:),POINTER ::dp
TYPE(dynadata),POINTER,INTENT(IN) :: head
CHARACTER(LEN=*),INTENT(IN) :: aname
!
TYPE(dynadata),POINTER :: iter=>NULL()
NULLIFY(dp)
iter => dyd_findnode(head,aname)
dp => iter%mreal8d4p
END SUBROUTINE
SUBROUTINE getdata_real8p5(head,aname,dp)
REAL(8),DIMENSION(:,:,:,:,:),POINTER ::dp
TYPE(dynadata),POINTER,INTENT(IN) :: head
CHARACTER(LEN=*),INTENT(IN) :: aname
!
TYPE(dynadata),POINTER :: iter=>NULL()
NULLIFY(dp)
iter => dyd_findnode(head,aname)
dp => iter%mreal8d5p
END SUBROUTINE
SUBROUTINE dyd_empty(this)
TYPE(dynadata),INTENT(INOUT) :: this
! local
LOGICAL fk
INTEGER ierr
this%mdim=-1; this%mrange=-1; this%metype=-1; this%mname=''
NULLIFY(this%mdynadata)
NULLIFY(this%mprev)
NULLIFY(this%mnext)
NULLIFY(this%mhead)
IF(ASSOCIATED(this%mint4d0p)) THEN
DEALLOCATE(this%mint4d0p,STAT=ierr); NULLIFY(this%mint4d0p)
ENDIF
IF(ASSOCIATED(this%mint4d1p)) THEN
DEALLOCATE(this%mint4d1p,STAT=ierr)
NULLIFY(this%mint4d1p)
ENDIF
IF(ASSOCIATED(this%mint4d2p)) THEN
DEALLOCATE(this%mint4d2p,STAT=ierr)
NULLIFY(this%mint4d2p)
ENDIF
IF(ASSOCIATED(this%mint4d3p)) THEN
DEALLOCATE(this%mint4d3p,STAT=ierr)
NULLIFY(this%mint4d3p)
ENDIF
IF(ASSOCIATED(this%mint4d4p)) THEN
DEALLOCATE(this%mint4d4p,STAT=ierr)
NULLIFY(this%mint4d4p)
ENDIF
IF(ASSOCIATED(this%mint4d5p)) THEN
DEALLOCATE(this%mint4d5p,STAT=ierr)
NULLIFY(this%mint4d5p)
ENDIF
IF(ASSOCIATED(this%mreal8d0p)) THEN
DEALLOCATE(this%mreal8d0p,STAT=ierr); NULLIFY(this%mreal8d0p)
ENDIF
IF(ASSOCIATED(this%mreal8d1p)) THEN
DEALLOCATE(this%mreal8d1p,STAT=ierr); NULLIFY(this%mreal8d1p)
ENDIF
IF(ASSOCIATED(this%mreal8d2p)) THEN
DEALLOCATE(this%mreal8d2p,STAT=ierr); NULLIFY(this%mreal8d2p)
ENDIF
IF(ASSOCIATED(this%mreal8d3p)) THEN
DEALLOCATE(this%mreal8d3p,STAT=ierr); NULLIFY(this%mreal8d3p)
ENDIF
IF(ASSOCIATED(this%mreal8d4p)) THEN
DEALLOCATE(this%mreal8d4p,STAT=ierr); NULLIFY(this%mreal8d4p)
ENDIF
IF(ASSOCIATED(this%mreal8d5p)) THEN
DEALLOCATE(this%mreal8d5p,STAT=ierr); NULLIFY(this%mreal8d5p)
ENDIF
END SUBROUTINE
SUBROUTINE dyd_new(head,elementtype, dataname, ierr, &
range11, range12,range21, range22,range31, range32, &
range41, range42,range51, range52)
TYPE(dynadata),POINTER,INTENT(INOUT) :: head
CHARACTER(LEN=*),INTENT(IN) :: elementtype, dataname
INTEGER,INTENT(IN),OPTIONAL :: range11, range12,range21, range22, &
range31, range32, range41, range42,range51, range52
INTEGER,INTENT(OUT) :: ierr
! local
INTEGER :: atype,adim
INTEGER :: arange(2,10)
TYPE(dynadata),POINTER :: pdyd
CHARACTER(LEN=32) :: tstr
ierr=0
IF(.NOT.dyd_isvalidname(head,dataname)) THEN
ierr=1
RETURN
ENDIF
! get type
tstr = TRIM(elementtype)
CALL upper_case(tstr)
IF(tstr=='INTEGER(4)') THEN; atype=int4
ELSEIF(tstr=='REAL(8)') THEN; atype=real8
ELSE
ierr=2
RETURN
ENDIF
! get
adim=0
IF( PRESENT(range11) .AND. PRESENT(range12) ) THEN
IF(range11>=range12) THEN
WRITE(*,*) 'Invalid range'; ierr=3
RETURN
ENDIF
adim=1; arange(1,1)=range11; arange(2,1)=range12
! DIM 2
IF( PRESENT(range21) .AND. PRESENT(range22) ) THEN
IF(range21>=range22) THEN
WRITE(*,*) 'Invalid range'; ierr=3
RETURN
ENDIF
adim=2; arange(1,2)=range21; arange(2,2)=range22
! DIM 3
IF( PRESENT(range31) .AND. PRESENT(range32) ) THEN
IF(range31>=range32) THEN
WRITE(*,*) 'Invalid range'; ierr=3
RETURN
ENDIF
adim=3; arange(1,3)=range31; arange(2,3)=range32
! DIM 4
IF( PRESENT(range41) .AND. PRESENT(range42) ) THEN
IF(range41>=range42) THEN
WRITE(*,*) 'Invalid range'; ierr=3
RETURN
ENDIF
adim=4; arange(1,4)=range41; arange(2,4)=range42
! DIM 5
IF( PRESENT(range51) .AND. PRESENT(range52) ) THEN
IF(range51>=range52) THEN
WRITE(*,*) 'Invalid range'; ierr=3
RETURN
ENDIF
adim=5; arange(1,5)=range51; arange(2,5)=range52
ENDIF
ENDIF
ENDIF
ENDIF
ENDIF
!WRITE(*,*) 'dyd_new ',atype,adim,arange(:,1:adim),dataname
CALL dyd_new_0(head,atype,adim,arange,dataname)
END SUBROUTINE
SUBROUTINE dyd_new_1(head,elementtype, arrayrange, dataname, ierr)
TYPE(dynadata),POINTER,INTENT(INOUT) :: head
CHARACTER(LEN=*),INTENT(IN) :: elementtype, arrayrange,dataname
INTEGER,INTENT(OUT) :: ierr
! local
INTEGER :: atype,offset,suboffset,i,adim
INTEGER :: arange(2,10), loc(20),com(10)
TYPE(dynadata),POINTER :: pdyd
CHARACTER(LEN=32) :: tstr
!WRITE(*,*) 'new1 ',dataname
ierr=0
IF(.NOT.dyd_isvalidname(head,dataname)) THEN
ierr=1
RETURN
ENDIF
! get type
tstr = TRIM(elementtype)
CALL upper_case(tstr)
IF(tstr=='INTEGER(4)') THEN; atype=int4
ELSEIF(tstr=='REAL(8)') THEN; atype=real8
ELSE
ierr=2
RETURN
ENDIF
! get
adim=0; offset=1
DO
suboffset = INDEX(arrayrange(offset:),':')
IF(suboffset==0) EXIT
adim = adim+1
IF(adim>5) STOP
offset = offset+suboffset
ENDDO
!WRITE(*,*) 'adim=',adim,' ',arrayrange
IF(adim>0) THEN
loc(1) = INDEX(arrayrange,'(')
loc(adim+1) = INDEX(arrayrange,')')
com(1) = INDEX(arrayrange,':')
IF(adim>1) THEN
DO i=2,adim
loc(i) = loc(i-1)+INDEX(arrayrange(loc(i-1)+1:),',')
!write(*,*) loc(i),' ',arrayrange(loc(i-1)+1:)
com(i) = loc(i)+INDEX(arrayrange(loc(i)+1:),':')
!write(*,*) com(i),' ',arrayrange(loc(i)+1:)
ENDDO
ENDIF
write(*,*) loc(1:adim+1)
write(*,*) com(1:adim)
DO i=1,adim
READ( arrayrange(loc(i)+1 : com(i)-1 ), *) arange(1,i)
READ( arrayrange(com(i)+1 : loc(i+1)-1), *) arange(2,i)
ENDDO
ENDIF
CALL dyd_new_0(head,atype,adim,arange,dataname)
!NULLIFY(pdyd)
!ALLOCATE(pdyd)
!CALL dyd_construct(pdyd, atype,adim,arange,dataname)
!CALL dyd_insert(pdyd,head)
END SUBROUTINE
SUBROUTINE dyd_new_0(head,atype,adim,arange,aname)
TYPE(dynadata),POINTER,INTENT(INOUT) :: head
INTEGER,INTENT(IN) :: atype,adim
INTEGER,INTENT(IN),OPTIONAL :: arange(2,adim)
CHARACTER(LEN=*),INTENT(IN) :: aname
! local
TYPE(dynadata),POINTER :: pdyd
NULLIFY(pdyd)
IF(dyd_isvalidname(head,aname)) THEN
ALLOCATE(pdyd)
CALL dyd_construct(pdyd, atype,adim,arange,aname)
CALL dyd_insert(pdyd,head)
ENDIF
END SUBROUTINE
SUBROUTINE dyd_construct(this,atype,adim,arange,aname)
TYPE(dynadata),INTENT(INOUT) :: this
INTEGER,INTENT(IN) :: atype,adim
INTEGER,INTENT(IN),OPTIONAL :: arange(2,adim)
CHARACTER(LEN=*),INTENT(IN) :: aname
!
INTEGER i
CALL dyd_empty(this)
IF(atype<=0 .OR. atype>2) THEN
WRITE(*,*) 'Invalid type: ',atype
STOP
ENDIF
IF(adim<0) THEN
WRITE(*,*) 'Invalid dim: ',adim
STOP
ENDIF
IF(adim>0) THEN; DO i=1,adim
IF(arange(1,i)>arange(2,i)) THEN
WRITE(*,*) 'Invalid range: ',arange(1:2,i), ' in DIM ',i
STOP
ENDIF
ENDDO; ENDIF
this%mname = aname
this%metype=atype
this%mdim=adim
IF(adim>0) this%mrange(1:2,1:adim) = arange(1:2,1:adim)
SELECT CASE(this%metype)
CASE (int4)
IF(this%mdim==0) THEN
ALLOCATE(this%mint4d0p)
ENDIF
IF(this%mdim==1) THEN
ALLOCATE(this%mint4d1p(this%mrange(1,1):this%mrange(2,1) ) )
ENDIF
IF(this%mdim==2) THEN
ALLOCATE(this%mint4d2p(this%mrange(1,1):this%mrange(2,1), &
this%mrange(1,2):this%mrange(2,2) ) )
ENDIF
IF(this%mdim==3) THEN
ALLOCATE(this%mint4d3p(this%mrange(1,1):this%mrange(2,1), &
this%mrange(1,2):this%mrange(2,2), &
this%mrange(1,3):this%mrange(2,3) ) )
ENDIF
IF(this%mdim==4) THEN
ALLOCATE(this%mint4d4p(this%mrange(1,1):this%mrange(2,1), &
this%mrange(1,2):this%mrange(2,2), &
this%mrange(1,3):this%mrange(2,3), &
this%mrange(1,4):this%mrange(2,4) ) )
ENDIF
IF(this%mdim==5) THEN
ALLOCATE(this%mint4d5p(this%mrange(1,1):this%mrange(2,1), &
this%mrange(1,2):this%mrange(2,2), &
this%mrange(1,3):this%mrange(2,3), &
this%mrange(1,4):this%mrange(2,4), &
this%mrange(1,5):this%mrange(2,5) ) )
ENDIF
CASE (real8)
IF(this%mdim==0) THEN
ALLOCATE(this%mreal8d0p)
ENDIF
IF(this%mdim==1) THEN
ALLOCATE(this%mreal8d1p(this%mrange(1,1):this%mrange(2,1) ) )
ENDIF
IF(this%mdim==2) THEN
ALLOCATE(this%mreal8d2p(this%mrange(1,1):this%mrange(2,1), &
this%mrange(1,2):this%mrange(2,2) ) )
ENDIF
IF(this%mdim==3) THEN
ALLOCATE(this%mreal8d3p(this%mrange(1,1):this%mrange(2,1), &
this%mrange(1,2):this%mrange(2,2), &
this%mrange(1,3):this%mrange(2,3) ) )
ENDIF
IF(this%mdim==4) THEN
ALLOCATE(this%mreal8d4p(this%mrange(1,1):this%mrange(2,1), &
this%mrange(1,2):this%mrange(2,2), &
this%mrange(1,3):this%mrange(2,3), &
this%mrange(1,4):this%mrange(2,4) ) )
ENDIF
IF(this%mdim==5) THEN
ALLOCATE(this%mreal8d5p(this%mrange(1,1):this%mrange(2,1), &
this%mrange(1,2):this%mrange(2,2), &
this%mrange(1,3):this%mrange(2,3), &
this%mrange(1,4):this%mrange(2,4), &
this%mrange(1,5):this%mrange(2,5) ) )
ENDIF
END SELECT
END SUBROUTINE
SUBROUTINE dyd_insert(this,head)
TYPE(dynadata),POINTER,INTENT(INOUT) :: this
TYPE(dynadata),POINTER,INTENT(INOUT) :: head
!
TYPE(dynadata),POINTER :: tail
IF( ASSOCIATED(head) ) THEN
CALL dyd_gettail(head,tail)
!WRITE(*,*) 'Tail name is ',TRIM(tail%mname)
this%mprev => tail
this%mhead => head
tail%mnext => this
NULLIFY(this%mnext)
ELSE
head => this
this%mprev => NULL()
this%mhead => this
this%mnext => NULL()
ENDIF
END SUBROUTINE
SUBROUTINE dyd_gettail(this,tail)
TYPE(dynadata),POINTER,INTENT(IN) :: this
TYPE(dynadata),POINTER,INTENT(OUT) :: tail
IF(ASSOCIATED(this)) THEN
tail => this
DO WHILE( ASSOCIATED(tail%mnext) )
tail => tail%mnext
END DO
ENDIF
END SUBROUTINE
FUNCTION dyd_isvalidname(head,aname)
TYPE(dynadata),POINTER,INTENT(IN) :: head
CHARACTER(LEN=*),INTENT(IN) :: aname
LOGICAL dyd_isvalidname
!
TYPE(dynadata),POINTER :: iter
dyd_isvalidname = .TRUE.
IF( ASSOCIATED(head) ) THEN
iter=>head
DO WHILE( ASSOCIATED(iter%mnext) )
!WRITE(*,*) TRIM(aname), ' ',TRIM(iter%mname)
!READ(*,*)
IF(TRIM(aname)==TRIM(iter%mname)) THEN
dyd_isvalidname = .FALSE.
EXIT
ENDIF
iter => iter%mnext
END DO
ENDIF
END FUNCTIOn
FUNCTION dyd_findnode(head,aname)
TYPE(dynadata),POINTER,INTENT(IN) :: head
CHARACTER(LEN=*),INTENT(IN) :: aname
TYPE(dynadata),POINTER :: dyd_findnode
!
TYPE(dynadata),POINTER :: iter=>NULL()
NULLIFY(dyd_findnode)
IF(.NOT.ASSOCIATED(head)) RETURN
iter=>head
loop: DO
IF(TRIM(aname)==TRIM(iter%mname)) THEN
dyd_findnode => iter ! 找到
EXIT loop
ENDIF
IF(.NOT.ASSOCIATED(iter%mnext)) THEN
EXIT loop ! 队列结束，没找到
ELSE
iter=>iter%mnext
ENDIF
END DO loop
END FUNCTION
! other
elemental subroutine upper_case(word)
! convert a word to lower case
character (len=*) , intent(in out) :: word
integer :: i,ic,nlen
nlen = len(word)
do i=1,nlen
ic = ichar(word(i:i))
if (ic >= 97 .and. ic < 122) word(i:i) = char(ic-32)
end do
end subroutine upper_case
END MODULE
测试代码如下
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
SUBROUTINE TEST_int4d0
USE xx
IMPLICIT NONE
TYPE(dynadata),POINTER :: list=>NULL(),iter=>NULL()
INTEGER(4),POINTER :: X
INTEGER ::ierr
!CALL dyd_new(head=list, atype=int4, adim=1,arange=irange, aname='X')
!CALL dyd_new(head=list, elementtype='integer(4)', arrayrange='(0:9)', &
! dataname='X', ierr=ierr)
CALL dyd_new(head=list, elementtype='integer(4)', dataname='X', ierr=ierr)
CALL getdata(list,'X', X)
iter => dyd_findnode(head=list,aname='X'); !X => iter%mint4d1p
WRITE(*,*) X
X = 5
WRITE(*,*)iter%mint4d0p
CALL dyd_empty(iter)
END SUBROUTINE
SUBROUTINE TEST_int4d1
USE xx
IMPLICIT NONE
TYPE(dynadata),POINTER :: list=>NULL(),iter=>NULL()
INTEGER(4),DIMENSION(:),POINTER :: X
INTEGER ::i, irange(2,2),ierr
irange(1,1) =0; irange(2,1) =9
!CALL dyd_new(head=list, atype=int4, adim=1,arange=irange, aname='X')
!CALL dyd_new(head=list, elementtype='integer(4)', arrayrange='(0:9)', &
! dataname='X', ierr=ierr)
CALL dyd_new(head=list, elementtype='integer(4)', dataname='X', ierr=ierr, &
range11=0,range12=9)
CALL getdata(list,'X', X)
iter => dyd_findnode(head=list,aname='X'); !X => iter%mint4d1p
WRITE(*,*) X
WRITE(*,*) SHAPE(X),LBOUND(X),UBOUND(X)
X = 5
DO i=0,9; X(i) = i; ENDDO
WRITE(*,*)iter%mint4d1p
CALL dyd_empty(iter)
END SUBROUTINE
SUBROUTINE TEST_int4d2
USE xx
IMPLICIT NONE
TYPE(dynadata),POINTER :: list=>NULL(),iter=>NULL()
INTEGER(4),DIMENSION(:,:),POINTER :: X
INTEGER ::i,j,k, irange(2,2),ierr
irange(1,1) =0; irange(2,1) =9
irange(1,2) =1; irange(2,2) =5
!CALL dyd_new(head=list, atype=int4, adim=1,arange=irange, aname='X')
!CALL dyd_new(head=list, elementtype='integer(4)', arrayrange='(0:9)', &
! dataname='X', ierr=ierr)
CALL dyd_new(head=list, elementtype='integer(4)', dataname='X', ierr=ierr, &
range11=0,range12=9,range21=1,range22=5)
CALL getdata(list,'X', X)
iter => dyd_findnode(head=list,aname='X'); !X => iter%mint4d1p
WRITE(*,*) X
WRITE(*,*) 'SHAPE(X)=',SHAPE(X),'LBOUND(X)=',LBOUND(X),'UBOUND(X)=',UBOUND(X)
k=1
DO j=1,5;DO i=0,9; X(i,j) = k; k=k+1; ENDDO;ENDDO
WRITE(*,*)iter%mint4d2p
CALL dyd_empty(iter)
END SUBROUTINE
SUBROUTINE TEST_int4d3
USE xx
IMPLICIT NONE
TYPE(dynadata),POINTER :: list=>NULL(),iter=>NULL()
INTEGER(4),DIMENSION(:,:,:),POINTER :: X
INTEGER ::i,j,k,l,m,n, irange(2,3),ierr
irange(1,1) =0; irange(2,1) =9
irange(1,2) =1; irange(2,2) =5
irange(1,3) =1; irange(2,3) =3
!CALL dyd_new(head=list, atype=int4, adim=3,arange=irange, aname='X')
CALL dyd_new(head=list, elementtype='integer(4)', arrayrange='(0:9,1:5,1:3)', &
dataname='X', ierr=ierr)
!CALL dyd_new(head=list, elementtype='integer(4)', dataname='X', ierr=ierr, &
! range11=0,range12=9,range21=1,range22=5,range31=1,range32=3)
CALL getdata(list,'X', X)
iter => dyd_findnode(head=list,aname='X'); !X => iter%mint4d1p
WRITE(*,*) X
WRITE(*,*) 'SHAPE(X)=',SHAPE(X),'LBOUND(X)=',LBOUND(X),'UBOUND(X)=',UBOUND(X)
k=1
DO l=1,3; DO j=1,5;DO i=0,9
X(i,j,l) = k; k=k+1
ENDDO;ENDDO;ENDDO
WRITE(*,*)iter%mint4d3p
CALL dyd_empty(iter)
END SUBROUTINE
SUBROUTINE TEST_int4d4
USE xx
IMPLICIT NONE
TYPE(dynadata),POINTER :: list=>NULL(),iter=>NULL()
INTEGER(4),DIMENSION(:,:,:,:),POINTER :: X
INTEGER ::i,j,k,l,m,n, irange(2,4),ierr
irange(1,1) =0; irange(2,1) =2
irange(1,2) =1; irange(2,2) =5
irange(1,3) =1; irange(2,3) =3
irange(1,4) =-1; irange(2,4) =1
CALL dyd_new(head=list, atype=int4, adim=4,arange=irange, aname='X')
!CALL dyd_new(head=list, elementtype='integer(4)', arrayrange='(0:9,1:5,1:3)', &
! dataname='X', ierr=ierr)
!CALL dyd_new(head=list, elementtype='integer(4)', dataname='X', ierr=ierr, &
! range11=0,range12=9,range21=1,range22=5,range31=1,range32=3)
CALL getdata(list,'X', X)
iter => dyd_findnode(head=list,aname='X'); !X => iter%mint4d1p
WRITE(*,*) X
WRITE(*,*) 'SHAPE(X)=',SHAPE(X),'LBOUND(X)=',LBOUND(X),'UBOUND(X)=',UBOUND(X)
k=1
DO m=irange(1,4), irange(2,4)
DO l=irange(1,3), irange(2,3)
DO j=irange(1,2), irange(2,2)
DO i=irange(1,1), irange(2,1)
X(i,j,l,m) = k; k=k+1
ENDDO;ENDDO;ENDDO;ENDDO
WRITE(*,*)iter%mint4d4p
CALL dyd_empty(iter)
END SUBROUTINE
SUBROUTINE TEST_int4d5
USE xx
IMPLICIT NONE
TYPE(dynadata),POINTER :: list=>NULL(),iter=>NULL()
INTEGER(4),DIMENSION(:,:,:,:,:),POINTER :: X
INTEGER ::i,j,k,l,m,n, irange(2,5),ierr
irange(1,1) =0; irange(2,1) =2
irange(1,2) =1; irange(2,2) =5
irange(1,3) =1; irange(2,3) =3
irange(1,4) =-1; irange(2,4) =1
irange(1,5) =5; irange(2,5) =8
!CALL dyd_new(head=list, atype=int4, adim=5,arange=irange, aname='X')
CALL dyd_new(head=list, elementtype='integer(4)', arrayrange='(0:2,1:5,1:3,-1:1,5:8)', &
dataname='X', ierr=ierr)
!CALL dyd_new(head=list, elementtype='integer(4)', dataname='X', ierr=ierr, &
! range11=0,range12=9,range21=1,range22=5,range31=1,range32=3)
CALL getdata(list,'X', X)
iter => dyd_findnode(head=list,aname='X'); !X => iter%mint4d1p
!WRITE(*,*) X
WRITE(*,*) 'SHAPE(X)=',SHAPE(X),'LBOUND(X)=',LBOUND(X),'UBOUND(X)=',UBOUND(X)
k=1
DO n=irange(1,5), irange(2,5)
DO m=irange(1,4), irange(2,4)
DO l=irange(1,3), irange(2,3)
DO j=irange(1,2), irange(2,2)
DO i=irange(1,1), irange(2,1)
X(i,j,l,m,n) = k; k=k+1
ENDDO;ENDDO;ENDDO;ENDDO;ENDDO
WRITE(*,*)iter%mint4d5p
CALL dyd_empty(iter)
END SUBROUTINE
SUBROUTINE TEST_real8d0
USE xx
IMPLICIT NONE
TYPE(dynadata),POINTER :: list=>NULL(),iter=>NULL()
REAL(8),POINTER :: X
INTEGER ::ierr
!CALL dyd_new(head=list, atype=real8, adim=1,arange=irange, aname='X')
!CALL dyd_new(head=list, elementtype='REAL(8)', arrayrange='(0:9)', &
! dataname='X', ierr=ierr)
CALL dyd_new(head=list, elementtype='REAL(8)', dataname='X', ierr=ierr)
CALL getdata(list,'X', X)
iter => dyd_findnode(head=list,aname='X'); !X => iter%mreal8d1p
!WRITE(*,*) X
X = 5
WRITE(*,*)iter%mreal8d0p
CALL dyd_empty(iter)
END SUBROUTINE
SUBROUTINE TEST_real8d1
USE xx
IMPLICIT NONE
TYPE(dynadata),POINTER :: list=>NULL(),iter=>NULL()
REAL(8),DIMENSION(:),POINTER :: X
INTEGER ::i, irange(2,2),ierr
irange(1,1) =0; irange(2,1) =9
!CALL dyd_new(head=list, atype=real8, adim=1,arange=irange, aname='X')
!CALL dyd_new(head=list, elementtype='REAL(8)', arrayrange='(0:9)', &
! dataname='X', ierr=ierr)
CALL dyd_new(head=list, elementtype='REAL(8)', dataname='X', ierr=ierr, &
range11=0,range12=9)
CALL getdata(list,'X', X)
iter => dyd_findnode(head=list,aname='X'); !X => iter%mreal8d1p
!WRITE(*,*) X
WRITE(*,*) SHAPE(X),LBOUND(X),UBOUND(X)
X = 5
DO i=0,9; X(i) = i; ENDDO
WRITE(*,*)iter%mreal8d1p
CALL dyd_empty(iter)
END SUBROUTINE
SUBROUTINE TEST_real8d2
USE xx
IMPLICIT NONE
TYPE(dynadata),POINTER :: list=>NULL(),iter=>NULL()
REAL(8),DIMENSION(:,:),POINTER :: X
INTEGER ::i,j,k, irange(2,2),ierr
irange(1,1) =0; irange(2,1) =9
irange(1,2) =1; irange(2,2) =5
!CALL dyd_new(head=list, atype=real8, adim=1,arange=irange, aname='X')
!CALL dyd_new(head=list, elementtype='REAL(8)', arrayrange='(0:9)', &
! dataname='X', ierr=ierr)
CALL dyd_new(head=list, elementtype='REAL(8)', dataname='X', ierr=ierr, &
range11=0,range12=9,range21=1,range22=5)
CALL getdata(list,'X', X)
iter => dyd_findnode(head=list,aname='X'); !X => iter%mreal8d1p
!WRITE(*,*) X
WRITE(*,*) 'SHAPE(X)=',SHAPE(X),'LBOUND(X)=',LBOUND(X),'UBOUND(X)=',UBOUND(X)
k=1
DO j=1,5;DO i=0,9; X(i,j) = k; k=k+1; ENDDO;ENDDO
WRITE(*,*)iter%mreal8d2p
CALL dyd_empty(iter)
END SUBROUTINE
SUBROUTINE TEST_real8d3
USE xx
IMPLICIT NONE
TYPE(dynadata),POINTER :: list=>NULL(),iter=>NULL()
REAL(8),DIMENSION(:,:,:),POINTER :: X
INTEGER ::i,j,k,l,m,n, irange(2,3),ierr
irange(1,1) =0; irange(2,1) =9
irange(1,2) =1; irange(2,2) =5
irange(1,3) =1; irange(2,3) =3
!CALL dyd_new(head=list, atype=real8, adim=3,arange=irange, aname='X')
CALL dyd_new(head=list, elementtype='REAL(8)', arrayrange='(0:9,1:5,1:3)', &
dataname='X', ierr=ierr)
!CALL dyd_new(head=list, elementtype='REAL(8)', dataname='X', ierr=ierr, &
! range11=0,range12=9,range21=1,range22=5,range31=1,range32=3)
CALL getdata(list,'X', X)
iter => dyd_findnode(head=list,aname='X'); !X => iter%mreal8d1p
!WRITE(*,*) X
WRITE(*,*) 'SHAPE(X)=',SHAPE(X),'LBOUND(X)=',LBOUND(X),'UBOUND(X)=',UBOUND(X)
k=1
DO l=1,3; DO j=1,5;DO i=0,9
X(i,j,l) = k; k=k+1
ENDDO;ENDDO;ENDDO
WRITE(*,*)iter%mreal8d3p
CALL dyd_empty(iter)
END SUBROUTINE
SUBROUTINE TEST_real8d4
USE xx
IMPLICIT NONE
TYPE(dynadata),POINTER :: list=>NULL(),iter=>NULL()
REAL(8),DIMENSION(:,:,:,:),POINTER :: X
INTEGER ::i,j,k,l,m,n, irange(2,4),ierr
irange(1,1) =0; irange(2,1) =2
irange(1,2) =1; irange(2,2) =5
irange(1,3) =1; irange(2,3) =3
irange(1,4) =-1; irange(2,4) =1
CALL dyd_new(head=list, atype=real8, adim=4,arange=irange, aname='X')
!CALL dyd_new(head=list, elementtype='REAL(8)', arrayrange='(0:9,1:5,1:3)', &
! dataname='X', ierr=ierr)
!CALL dyd_new(head=list, elementtype='REAL(8)', dataname='X', ierr=ierr, &
! range11=0,range12=9,range21=1,range22=5,range31=1,range32=3)
CALL getdata(list,'X', X)
iter => dyd_findnode(head=list,aname='X'); !X => iter%mreal8d1p
!WRITE(*,*) X
WRITE(*,*) 'SHAPE(X)=',SHAPE(X),'LBOUND(X)=',LBOUND(X),'UBOUND(X)=',UBOUND(X)
k=1
DO m=irange(1,4), irange(2,4)
DO l=irange(1,3), irange(2,3)
DO j=irange(1,2), irange(2,2)
DO i=irange(1,1), irange(2,1)
X(i,j,l,m) = k; k=k+1
ENDDO;ENDDO;ENDDO;ENDDO
WRITE(*,*)iter%mreal8d4p
CALL dyd_empty(iter)
END SUBROUTINE
SUBROUTINE TEST_real8d5
USE xx
IMPLICIT NONE
TYPE(dynadata),POINTER :: list=>NULL(),iter=>NULL()
REAL(8),DIMENSION(:,:,:,:,:),POINTER :: X
INTEGER ::i,j,k,l,m,n, irange(2,5),ierr
irange(1,1) =0; irange(2,1) =2
irange(1,2) =1; irange(2,2) =5
irange(1,3) =1; irange(2,3) =3
irange(1,4) =-1; irange(2,4) =1
irange(1,5) =5; irange(2,5) =8
!CALL dyd_new(head=list, atype=real8, adim=5,arange=irange, aname='X')
CALL dyd_new(head=list, elementtype='REAL(8)', arrayrange='(0:2,1:5,1:3,-1:1,5:8)', &
dataname='X', ierr=ierr)
!CALL dyd_new(head=list, elementtype='REAL(8)', dataname='X', ierr=ierr, &
! range11=0,range12=9,range21=1,range22=5,range31=1,range32=3)
CALL getdata(list,'X', X)
iter => dyd_findnode(head=list,aname='X'); !X => iter%mreal8d1p
!WRITE(*,*) X
WRITE(*,*) 'SHAPE(X)=',SHAPE(X),'LBOUND(X)=',LBOUND(X),'UBOUND(X)=',UBOUND(X)
k=1
DO n=irange(1,5), irange(2,5)
DO m=irange(1,4), irange(2,4)
DO l=irange(1,3), irange(2,3)
DO j=irange(1,2), irange(2,2)
DO i=irange(1,1), irange(2,1)
X(i,j,l,m,n) = k; k=k+1
ENDDO;ENDDO;ENDDO;ENDDO;ENDDO
WRITE(*,*)iter%mreal8d5p
CALL dyd_empty(iter)
END SUBROUTINE
PROGRAM TEST
IMPLICIT NONE
! CALL TEST_int4d0
! CALL TEST_int4d1
! CALL TEST_int4d2
! CALL TEST_int4d3
! CALL TEST_int4d4
! CALL TEST_int4d5
CALL TEST_real8d0
CALL TEST_real8d1
CALL TEST_real8d2
CALL TEST_real8d3
CALL TEST_real8d4
CALL TEST_real8d5
END
FEM相关
SUPG中的稳定项是二阶偏导，人工耗散系数与ah成正比，a是特征速度，h是网格尺度。
JST中心差分的稳定项是四阶偏导。
因此我以为JST中心差分的稳定项的阶段误差更低。
但是仔细想想，JST人工耗散系数是谱半径成正比，考虑到除单元体积，则相当于aA/Vol，其中A/Vol就是网格尺度h。因此，SUPG的人工耗散系数和Jameson中心差分的是一样。只是梯度阶数不一样，一个是二阶，一个是四阶。但是两者带来截断误差都是O(h)。
扫了一眼开发者的博士论文论文，LM的网格存储没有按照分区进行完全分布式存储，而是每个进程都保留全部网格，这样就限制了并行的规模。作者当时进行的最大规模并行计算是400个进程。相比而言，DII的并行能力要成熟许多，测试过数千进程。如果想做FEM开发，似乎DII还是更好的选择。