incomp-flame-post/code/code_gen/post.py
ignis 1538773c83
All checks were successful
Deploy Documentation / build-and-deploy (push) Successful in 55s
docs: extensively document code generator pipeline stages and sympy optimization in post.py
2026-06-02 04:03:19 +00:00

1742 lines
68 KiB
Python

"""
DNS Post-Processing Code Generator Core (post.py)
==================================================
이 모듈은 난류 및 연소 DNS(Direct Numerical Simulation) 데이터의 후처리를 위한 고성능 Fortran 코드를 생성하는 컴파일러의 코어입니다.
사용자가 정의한 DSL(Domain Specific Language) 입력 식을 읽어 파싱한 후, 다음과 같은 4단계 최적화 컴파일 과정을 거쳐 극도로 최적화된 3차원 루프 Fortran 모듈 코드를 자동 생성합니다.
[컴파일러 파이프라인 4단계 개요]
1. Stage 1 (AST 수집 및 변수 정의):
- Lark 파서가 생성한 AST(Abstract Syntax Tree)를 순회하며 기본 입력 변수(Primary), 계산이 필요한 대입식(Derived),
그리고 통계 물리량 계산을 위한 평균화 변수(Averaged)를 추출하고 필드 메타데이터 객체를 구성합니다.
2. Stage 2 (수치 미분 및 변동량 확장):
- 수식 내에 수치 미분자(ddx, d2dy 등)나 변동량(fluctuation, u')이 존재하면, 이를 물리적으로 차분 연산할 중간
미분 필드(DerivedField) 및 변동량 필드(FluctuationField)로 자동 변환하고 변수 테이블에 등록하여 확장합니다.
3. Stage 3 (의존성 분석 및 위상 정렬):
- 필드 간의 선후 연산 관계를 분석하여 유향 의존성 그래프(Directed Dependency Graph)를 생성합니다.
- 난류 통계의 특성상 평균 연산을 기준으로 "평균치 계산 전의 루프(Pass 1)""평균치를 구한 후 변동량을 계산하는 루프(Pass 2)"
전체 연산 블록을 논리적으로 분할하고, 각각의 블록 내에서 올바른 순서로 계산되도록 위상 정렬(Topological Sort)을 수행합니다.
4. Stage 4 (수식 기호 최적화, Liveness 분석 및 Buffer Array Pooling):
- SymPy 기호 수학 라이브러리를 이용하여 복잡한 3차원 수식을 대수적으로 간소화하고, 공통 부분 식 제거(CSE)를 적용하여 연산 비용(Flops)을 최적화합니다.
- 격자 데이터가 거대하므로 모든 변수에 개별 3D 배열을 할당하면 메모리가 고갈됩니다. 이를 방지하기 위해 변수들의
생명 주기(Liveness Window)를 수학적으로 추적하고, 동적 메모리 풀을 구축하여 동시에 활성화되지 않는 임시 변수들이
공통의 제한된 버퍼 배열(xyzbuffer0, xyzbuffer1, ...)을 나누어 사용(Array Pooling)하도록 할당하여 메모리 사용량을 최소화합니다.
"""
from lark import Lark, Visitor, Transformer, v_args, Token
import warnings
from jinja2 import Template
import sympy
from sympy.printing.fortran import FCodePrinter
@v_args(inline=True)
class LarkToSympy(Transformer):
"""Lark AST의 수학적 노드들을 SymPy 기호 수식 객체로 변환하는 Transformer 클래스입니다.
이를 통해 DSL 수식 트리를 SymPy 기호 객체로 매핑하고, SymPy의 강력한 수식 간소화(Simplify)
및 공통 부분식 제거(CSE) 최적화 도구를 파이프라인 하부에서 활용할 수 있게 합니다.
"""
def __init__(self, fdict):
"""Lark-to-SymPy Transformer를 초기화합니다.
Args:
fdict (dict): 변수명에서 물리 필드 정의 객체(FieldBase 자식 클래스들)로 매핑되는 딕셔너리.
"""
self.fdict = fdict
def number(self, numeral):
# 상수를 SymPy Float 기호 객체로 변환
return sympy.Float(float(numeral))
def env(self, name):
# 환경 변수($ 기호로 시작)를 SymPy Symbol로 매핑
return sympy.Symbol(name.value)
def paren(self, val):
# 괄호는 수식 우선순위 보존 후 내부 식 반환
return val
def var(self, name):
# 일반 변수명을 SymPy Symbol로 매핑
return sympy.Symbol(name.value)
def fluc(self, name):
# 변동량 u' 기호를 SymPy가 해석할 수 있도록 u__prime 기호명으로 매핑
return sympy.Symbol(name.value + "__prime")
def dnx(self, partial, b):
# ddx(u)와 같은 수치 미분항을 ddx_u 형태의 하나의 고유 Symbol로 묶어서 취급
signature = f"{partial.data}_{b.value}"
return sympy.Symbol(signature)
def icall(self, op, val):
# sqr, pow3 등의 인라인 전용 함수들을 SymPy 지수 표현식으로 다이렉트 변환
if op.data == "sqr":
return val**2
elif op.data == "pow3":
return val**3
return val
def fcall(self, *args):
# 내장 수학 함수(sqrt, exp, log, abs, rxn_rate) 또는 사용자 정의 UDF들을 SymPy 함수 노드로 매핑
a = args[0]
func_name = a.value if hasattr(a, 'value') else str(a)
if func_name == "sqrt":
return sympy.sqrt(args[1])
elif func_name == "exp":
return sympy.exp(args[1])
elif func_name == "log":
return sympy.log(args[1])
elif func_name == "abs":
return sympy.Abs(args[1])
elif func_name == "rxn_rate":
return sympy.Function("rxn_rate")(args[1])
elif func_name == "udf":
return sympy.Function(a.value)(*args[1:])
return sympy.Function(func_name)(*args[1:])
def neg(self, val):
return -val
def add(self, a, b):
return a + b
def sub(self, a, b):
return a - b
def mul(self, a, b):
return a * b
def div(self, a, b):
return a / b
def udf(self, a):
return a.value
log = lambda self: "log"
exp = lambda self: "exp"
sqrt = lambda self: "sqrt"
abs = lambda self: "abs"
rxn_rate = lambda self: "rxn_rate"
class ArrayFCodePrinter(FCodePrinter):
"""SymPy 식을 Fortran 코드로 인쇄할 때, 변수들을 적절한 다차원 격자점 인덱스 배열(예: var(i,j,k))
또는 평균 배열(예: var(i)) 형식으로 포맷팅하여 변환 출력해 주는 특수 Printer 클래스입니다.
"""
def __init__(self, settings=None, array_symbols=None, avg_symbols=None):
settings = settings or {}
settings.setdefault('source_format', 'free') # 기본 프리포맷 Fortran 95 준수
settings.setdefault('standard', 95)
super().__init__(settings)
self.array_symbols = array_symbols or {} # 3차원 필드 배열의 실제 물리 버퍼 명칭 매핑
self.avg_symbols = avg_symbols or {} # 평균 필드 배열의 명칭 매핑
def _print_Float(self, expr):
# 배정밀도 실수 상수에 d0 접미사를 명시하여 컴파일러 정밀도 유실 방지
val = str(expr)
if 'e' in val or 'E' in val:
return val.replace('e', 'd').replace('E', 'd')
if '.' not in val:
return val + ".0d0"
return val + "d0"
def _print_Symbol(self, expr):
name = expr.name
# 1. 3D 공간 격자점 버퍼 할당명에 매핑되어 있는 경우 (i,j,k) 차분 격자점 인덱스를 부착
if name in self.array_symbols:
return f"{self.array_symbols[name]}(i,j,k)"
# 2. X방향(i) 1차원 평균화 변수인 경우 (i) 인덱스 부착
if name in self.avg_symbols:
return f"{self.avg_symbols[name]}(i)"
if name.startswith("avg_") or "_avg_" in name or name.endswith("_avg"):
return f"{name}(i)"
# 3. 변동량 변수(u__prime)인 경우, 원본 격자점 값과 1차원 평균값의 편차식인 '(u(i,j,k) - avg_u(i))'로 자동 인라인 전개
if name.endswith("__prime"):
base = name[:-7]
arr = self.array_symbols.get(base, base)
return f"({arr}(i,j,k) - {{0}}avg_{base}(i))"
return name
def _print_Function(self, expr):
# SymPy 내장 함수가 아닌 사용자 정의 함수(rxn_rate, udf 등)의 정상적인 Fortran 출력을 위한 예외 폴백 처리
try:
return super()._print_Function(expr)
except Exception:
args = ", ".join(self.doprint(arg) for arg in expr.args)
return f"{expr.func.__name__}({args})"
class SympyOptimizer:
"""컴파일러 내부에서 SymPy 수학 기호 연산 엔진을 제어하는 최적화 매니저 클래스입니다.
1. 수식 확장 대입(Substitution): 사용자가 정의한 물리 변수 중 디스크로 내보내지(export) 않는
임시 변수들의 수식을 참조하는 하부 수식에 재귀적으로 대입하여 거대한 하나의 식으로 확장합니다.
2. 대수적 간소화(Simplify & Cancel): SymPy를 사용하여 분모 분자 소거 및 삼각/지수 대수 간소화를 적용합니다.
3. 공통 부분식 제거(CSE - Common Subexpression Elimination): 복잡하게 얽힌 다차원 미분 수식 속에서
중복 계산되는 서브 식(예: 공통 항곱)을 탐지하여, 루프 외/내부에서 임시 스칼라 변수(x0, x1 등)에 선언한 후
루프 내 연산에서 한 번만 대입하여 연산하도록 루프 코드를 개조하여 FLOPS와 캐시 로드를 획기적으로 낮춥니다.
"""
_instance = None
@classmethod
def get_instance(cls, fdict):
if cls._instance is None or cls._instance.fdict is not fdict:
cls._instance = cls(fdict)
return cls._instance
def __init__(self, fdict):
self.fdict = fdict
self.sympy_cache = {} # 중복 연산 계산 방지를 위한 SymPy 식 캐시
self.exported_fields = set(
name for name, f in fdict.items()
if hasattr(f, 'attr') and f.attr.get('export')
) # 디스크로 파일 출력이 설정된 최종 목표 변수 집합
self.averaged_targets = set()
def set_averaged(self, averaged_dict):
# 평균화 대상이 되는 변수 집합을 수집
self.averaged_targets = {a.target for a in averaged_dict.values()}
def get_sympy_expr(self, name):
"""특정 변수 명칭에 대해 기호적으로 완전히 확장 대입된 SymPy Expression 객체를 빌드하여 캐싱합니다.
기본 입력 변수(Primary), 미분 변수(Derived), 평균 변수(Averaged), 변동량 변수(Fluctuation)는
치환하지 않고 독립적인 리프 기호(Symbol)로 둡니다. 그 외 중간에서만 활용되는 일반 임시 Derived 변수들은
의존하는 원본 수식들을 탐색하여 재귀적으로 인라인 대체(expand) 시킴으로써, SymPy 기호 엔진이
루프 수식 전체의 논리적 최적화를 전역적으로 수행할 수 있도록 합니다.
"""
if name in self.sympy_cache:
return self.sympy_cache[name]
field = self.fdict[name]
if hasattr(field, 'prime') and field.prime:
expr = sympy.Symbol(name)
self.sympy_cache[name] = expr
return expr
if hasattr(field, 'op'): # DerivedField (ddx 등 미분 결과)
expr = sympy.Symbol(name)
self.sympy_cache[name] = expr
return expr
if hasattr(field, 'weighted'): # AveragedField
expr = sympy.Symbol(name)
self.sympy_cache[name] = expr
return expr
if hasattr(field, 'field') and hasattr(field, 'w'): # FluctuationField
expr = sympy.Symbol(name)
self.sympy_cache[name] = expr
return expr
transformer = LarkToSympy(self.fdict)
expr = transformer.transform(field.exp)
# 재귀적으로 의존 중인 캐싱되지 않은 중간 계산 임시 변수들을 수식 본문으로 확장 치환(Substitution)
expanded_expr = expr
changed = True
while changed:
changed = False
free_syms = list(expanded_expr.free_symbols)
sub_dict = {}
for sym in free_syms:
sym_name = sym.name
if sym_name in self.fdict:
f = self.fdict[sym_name]
is_derived_field = hasattr(f, 'op')
is_averaged_field = hasattr(f, 'weighted')
is_primary_field = hasattr(f, 'prime') and f.prime
is_exported = sym_name in self.exported_fields
is_averaged_target = sym_name in self.averaged_targets
# 최종 출력 대상이거나 물리적 경계 처리가 필요한 핵심 변수가 아니라면,
# 기 기하학적 수식을 최적화하기 위해 현재 식 안으로 대입해 버립니다.
if not (is_derived_field or is_averaged_field or is_primary_field or is_exported or is_averaged_target):
sub_dict[sym] = self.get_sympy_expr(sym_name)
changed = True
if sub_dict:
expanded_expr = expanded_expr.subs(sub_dict)
self.sympy_cache[name] = expanded_expr
return expanded_expr
def calculate_flops_and_heavy(self, expr):
"""수식에서 수행되는 부동 소수점 연산(FLOPS) 및 계산 부하가 큰 연산(Division, Sqrt, Exp 등)의 개수를 평가합니다.
이를 통해 최적화 리포트에 연산 비용 증감률을 리포팅합니다.
"""
flops = 0
heavy = 0
for node in sympy.preorder_traversal(expr):
if isinstance(node, sympy.Add):
flops += len(node.args) - 1
elif isinstance(node, sympy.Mul):
flops += len(node.args) - 1
elif isinstance(node, sympy.Pow):
base, exp = node.args
if exp == 0.5 or exp == -0.5:
flops += 10 # 제곱근 및 역제곱근은 대략 10 flops 가중치 부여
heavy += 1
elif exp == -1:
flops += 4 # 역수 나눗셈은 4 flops 가중치 부여
heavy += 1
elif isinstance(exp, sympy.Integer):
val = abs(int(exp))
if val > 1:
flops += val - 1
else:
flops += 10
heavy += 1
elif isinstance(node, (sympy.Derivative, sympy.Function)):
name = node.func.__name__
if name == 'sqrt':
flops += 10
heavy += 1
elif name in ('exp', 'log', 'sin', 'cos', 'tan', 'rxn_rate'):
flops += 10 # 특수 수학 함수 연산은 10 flops 가중치 및 Heavy 연산자로 평가
heavy += 1
elif name == 'Abs':
flops += 1
else:
flops += 10
heavy += 1
return flops, heavy
def count_3d_loads(self, expr, three_d_arrays):
"""수식 내부에서 참조하는 3차원 격자점 배열의 총 메모리 로드 횟수를 카운트합니다.
메모리 대역폭 한계(Memory-bound)에 부딪히는 HPC 연산의 병목을 평가하는 데 중요한 지표입니다.
"""
count = 0
for node in sympy.preorder_traversal(expr):
if isinstance(node, sympy.Symbol) and node.name in three_d_arrays:
count += 1
return count
def optimize_field(self, name, alloc=None):
"""개별 물리 필드 수식을 SymPy를 사용해 최적화하고 공통 부분식(CSE) 코드를 추출합니다.
Args:
name (str): 최적화할 필드 변수명.
alloc (dict): 변수명에서 버퍼 배열명(xyzbufferN)으로 매핑되는 테이블.
Returns:
rhs (str): 최적화가 완료된 최종 Fortran 우변(Right-Hand Side) 코드 문자열.
cse_decls (list): CSE 추출 결과 생성된 임시 실수 변수 선언 리스트 (예: real(real64) :: x0).
cse_assigns (list): 루프 내부 대입에 사용될 임시 변수의 연산식 리스트.
"""
expr = self.get_sympy_expr(name)
three_d_arrays = {
k for k, v in self.fdict.items()
if hasattr(v, 'dim') and v.dim == ':,:,:'
}
# 최적화 전 비용 수집
before_flops, before_heavy = self.calculate_flops_and_heavy(expr)
before_loads = self.count_3d_loads(expr, three_d_arrays)
# 1차 대수적 간소화 및 약분 소거
simplified_expr = sympy.simplify(expr)
simplified_expr = sympy.cancel(simplified_expr)
# 기호명을 실제 Fortran 격자 배열명(xyzbufferN)으로 매핑하기 위한 준비
array_symbols = {}
for k, v in self.fdict.items():
if hasattr(v, 'array') and v.array:
array_symbols[k] = v.array
elif alloc and k in alloc:
array_symbols[k] = alloc[k]
else:
array_symbols[k] = k
printer = ArrayFCodePrinter(array_symbols=array_symbols)
# 2차: 공통 부분식 제거(CSE) 알고리즘 가동
# 복잡하게 중복 사용되는 수식 구조를 서브 트리로 분리하여 임시 기호(x0, x1 등)와 최종 축소 식으로 분할
replacements, reduced_exprs = sympy.cse(simplified_expr)
reduced_expr = reduced_exprs[0]
# 최적화 후 비용 계산
after_flops = 0
after_heavy = 0
after_loads = 0
for temp_var, temp_expr in replacements:
f_val, h_val = self.calculate_flops_and_heavy(temp_expr)
after_flops += f_val
after_heavy += h_val
after_loads += self.count_3d_loads(temp_expr, three_d_arrays)
f_val, h_val = self.calculate_flops_and_heavy(reduced_expr)
after_flops += f_val
after_heavy += h_val
after_loads += self.count_3d_loads(reduced_expr, three_d_arrays)
import sys
def pct_str(before, after):
if before == 0:
return "0.0%" if after == 0 else "+inf%"
diff = after - before
pct = (diff / before) * 100
return f"{pct:+.1f}%"
flops_pct = pct_str(before_flops, after_flops)
heavy_pct = pct_str(before_heavy, after_heavy)
loads_pct = pct_str(before_loads, after_loads)
# 성능 개선 등급 진단
if after_flops < before_flops * 0.5 or after_loads < before_loads * 0.5:
est_speedup = "Highly significant"
elif after_flops < before_flops or after_loads < before_loads:
est_speedup = "Moderate"
else:
est_speedup = "Minimal / Already optimal"
# 컴파일 중 최적화 상세 리포트를 stderr로 콘솔 화면에 로깅
sys.stderr.write(f"\n[SymPy Optimizer Report: {name}]\n")
sys.stderr.write(f"- Floating Point Ops : {before_flops} -> {after_flops} ({flops_pct})\n")
sys.stderr.write(f"- Heavy Ops (Div/Sqrt): {before_heavy} -> {after_heavy} ({heavy_pct})\n")
sys.stderr.write(f"- 3D Array Mem Reads : {before_loads} -> {after_loads} ({loads_pct})\n")
sys.stderr.write(f"=> Estimated Speedup in loop: {est_speedup}\n\n")
cse_decls = []
cse_assigns = []
# CSE 추출 결과를 Fortran 코드로 포매팅
if replacements:
for temp_var, temp_expr in replacements:
# 임시 변수 로컬 스칼라 선언 구문 작성
cse_decls.append(f"real(real64) :: {temp_var}")
# 루프 내에서 먼저 연산 후 할당하는 구문 작성
cse_assigns.append(f"{temp_var} = {printer.doprint(temp_expr)}")
rhs = printer.doprint(reduced_expr)
return rhs, cse_decls, cse_assigns
class CollectDefinitions(Visitor):
"""Lark AST의 전체 노드를 스캔하며 정의문(변수 리스트, 대입식, 평균화 조건)을 수집하여
컴파일러의 빌딩 블록 객체들로 구축하는 Visitor 클래스입니다 (Stage 1 핵심 작동부).
"""
def __init__ (self, primary, derived, averaged):
self.primary = primary # 기본 입력 물리 필드명 집합 (격자 파일 등에서 읽음)
self.derived = derived # 유도/계산식 필드 딕셔너리
self.averaged = averaged # 평균 연산 조건 딕셔너리
def varlist(self, tree):
# 괄호 안 [u, v, w] 같은 형식으로 표기된 기본(Primary) 변수 리스트 파싱
for v in tree.children:
self.primary.add(v.value)
# 의존성 테이블에 기본 변수 객체 등록
self.derived[v.value] = PrimaryField(v.value, self.derived)
def assign_var (self, tree):
# 등호(=) 연산자를 통한 변수 계산식 정의 노드 파싱
# 속성값(예: (export=true, latex="..."))이 존재하면 추출
if len(tree.children) > 2:
lval, lattr, rval = tree.children
else:
lval, rval = tree.children
lattr = None
attr_dict = {}
if lattr is not None:
for t in lattr.children:
k, v = t.children
attr_dict[k.value] = v.value
if lval.value in self.derived:
raise ValueError("duplicate definition of " + lval)
# Field 객체를 생성하여 derived 계산 테이블에 등록
self.derived[lval.value] = Field(lval.value, attr_dict, rval, self.derived)
def assign_avg_var (self, tree):
# avg 또는 avg$w {u, v, ...} 와 같은 평균화 지시어 파싱 및 등록
w = tree.children[0]
targets = tree.children[1:]
if (not w.children) or (w.children[0] is None):
self.averaged[None] = set([x.value for x in targets])
else:
self.averaged[w.children[0].value] = set([x.value for x in targets])
class ExpInspector(Visitor):
"""수식 AST 노드들을 순회하며 해당 식 내부의 변동량(Fluctuation) 포함 여부,
참조 중인 다른 물리 변수 종속성 목록(dep), 그리고 수치 미분항(deriv) 목록을 추출하는 Visitor 클래스입니다.
"""
def __init__(self):
self.fluctuation = False
self.dep = set([])
self.deriv = set([])
@classmethod
def inspect(cls, tree):
self = cls()
return self(tree)
def __call__(self, tree):
self.visit(tree)
return self.fluctuation, self.dep, self.deriv
def fluc(self, tree):
# 식 내부에 변동량 항이 존재함을 기록하고 종속성에 등록
self.fluctuation = True
self.dep.add(tree.children[0].value)
def var(self, tree):
self.dep.add(tree.children[0].value)
def dnx (self, tree):
# ddx(u) 와 같이 수치 미분 노드를 만나면 ddx_u 형태의 중간 미분 변수를 종속성 리스트에 기록
op, v = tree.children
deriv = "{}_{}".format(op.data, v.value)
self.dep.add(deriv)
self.deriv.add((op.data, v.value))
@v_args(inline=True) # 메서드의 서명 형식을 자식 노드 매개변수 나열형으로 변환
class ExpToLatex(Transformer):
"""DSL 수학 표현식 AST를 문서를 검토하기 편하게 LaTeX 문법 수식으로 번역하는 번역기입니다."""
def __init__(self, fdict):
self.fdict = fdict
def arithmatic_rooted(self, name):
try:
exproot = self.fdict[name].exp.data
except AttributeError:
exproot = "something_11fasq2afa3rfzsaerqw23"
return ((exproot == "add") or (exproot == "sub") or
(exproot == "mul") or (exproot == "div"))
def parenthise(self, name):
try:
latex = self.fdict[name].latex
latex_given = self.fdict[name].latex_given
except KeyError:
warnings.warn(name + " is not found")
latex = r"\mathrm{{{}}}".format(name)
latex_given = None
if self.arithmatic_rooted(name) and (latex_given is None):
latex = "(" + latex + ")"
return latex
def number(self, numeral):
return numeral
def env(self, name):
return r"\mathrm{{{}}}".format(name.value)
def paren(self, name):
return "({})".format(str(name))
def var(self, name):
return self.parenthise(name.value)
def fluc(self, name):
# LaTeX 수식 표기용 u'' 프라임 기호 렌더링
return self.parenthise(name.value)+ "''"
def dnx (self, partial, b):
fmt = r"\partial_{{{}}}"
coord = partial.data[-1]
op = fmt.format(coord + coord if len(partial.data) > 3 else coord)
signature = "{}_{}".format(partial.data, b.value)
try:
eq = self.fdict[signature].latex
except KeyError:
eq = op + self.parenthise(b.value)
warnings.warn(signature + " is not found: " + eq)
return eq
def icall (self, a, b):
if a.data == "sqr":
fcode = "({0})^2".format(b)
elif a.data == "pow3":
fcode = "({0})^3".format(b)
else:
fcode = "({0})".format(b)
return fcode
def fcall (self, *args):
a = args[0]
b = ", ".join(args[1:])
if a == 'sqrt':
fcode = r"\sqrt{{{}}}".format(b)
return fcode
elif a == 'abs':
fcode = r"\left| {} \right|".format(b)
return fcode
elif a.startswith("\\"):
fcode = r"{}{{({})}}".format(a, b)
return fcode
else:
fcode = r"\mathrm{{{}}}({})".format(a, b)
return fcode
def neg(self, b):
fcode = "(-{})".format(b)
return fcode
def add(self, a, b):
fcode = "{} + {}".format(a, b)
return fcode
def sub(self, a, b):
fcode = "{} - {}".format(a, b)
return fcode
def mul(self, a, b):
fcode = "{} {}".format(a, b)
return fcode
def div(self, a, b):
fcode = "{} / {}".format(a, b)
return fcode
log = lambda self : "\log"
exp = lambda self : "\exp"
sqrt = lambda self : "sqrt"
abs = lambda self : "abs"
rxn_rate = lambda self : "\omega"
udf = lambda self, a : a.value
@v_args(inline=True)
class ExpToCode(Transformer):
"""DSL 수학 수식 트리를 Fortran의 3차원 인덱스 코드 형태의 수식 문자열로 다이렉트 매핑하는 번역기입니다.
SymPy 최적화를 사용하지 않는 대체 연산 경로 등에서 보조적으로 사용됩니다.
"""
def __init__(self, fdict):
self.fdict = fdict
def number(self, numeral):
return str(float(numeral))
def env(self, name):
return name.value
def paren(self, name):
return "({})".format(str(name))
def var(self, name):
try:
arrname = self.fdict[name.value].array
except KeyError:
arrname = name.value
return arrname + "(i,j,k)"
def fluc(self, name):
try:
arrname = self.fdict[name.value].array
except KeyError:
arrname = name.value
fmt = "({0}(i,j,k) - {{0}}avg_{1}(i))"
return fmt.format(arrname, name.value)
def dnx (self, partial, b):
signature = "{}_{}".format(partial.data, b.value)
try:
arrname = self.fdict[signature].array
except KeyError:
arrname = signature
return arrname + "(i,j,k)"
def icall (self, a, b):
if a.data == "sqr":
fcode = "(({0})*({0}))".format(b)
elif a.data == "pow3":
fcode = "(({0})*({0})*({0}))".format(b)
else:
fcode = "({0})".format(b)
return fcode
def fcall (self, *args):
a = args[0]
b = ", ".join(args[1:])
fcode = "( {} ( {} ) )".format(a, b)
return fcode
def neg(self, b):
fcode = "( - {} )".format(b)
return fcode
def add(self, a, b):
fcode = "( {} + {} )".format(a, b)
return fcode
def sub(self, a, b):
fcode = "( {} - {} )".format(a, b)
return fcode
def mul(self, a, b):
fcode = "( {} * {} )".format(a, b)
return fcode
def div(self, a, b):
fcode = "( {} / {} )".format(a, b)
return fcode
log = lambda self : "log"
exp = lambda self : "exp"
sqrt = lambda self : "sqrt"
abs = lambda self : "abs"
rxn_rate = lambda self : "rxn_rate"
udf = lambda self, a : a.value
def make_allocate(name, shape, init_zero=True):
"""Fortran 배열을 동적 할당하고 예외 발생 시 프로세스를 안전하게 폭파시키는 할당 코드를 작성해 줍니다."""
alloc_str = f"allocate({name}({shape}), stat=ierr)\n"
alloc_str += f"if (ierr /= 0) then\n"
alloc_str += f" write(0,*) 'Error: allocation of {name} failed on process', myid\n"
alloc_str += f" call MPI_ABORT(MPI_COMM_TASK, 1, mpi_err)\n"
alloc_str += f"end if"
if init_zero:
alloc_str += f"\n{name} = 0."
return alloc_str
class FieldBase (object):
"""모든 물리 필드 객체의 최상위 기본 클래스로서 공통 데이터 구조 및 선언/할당 문자열 작성을 담당합니다."""
def __init__ (self, name, fdict):
self.name = name
self.array = name # 실제 메모리 버퍼 상의 지칭 배열 명칭
self.dep = set([]) # 종속된(의존 중인) 하위 필드 리스트
self.fluc = False # 변동량 편차 계산 대상인지 여부
self.prime = False # 기본 데이터 파일로부터 읽어 들이는 원본 필드인지 여부
self.fdict = fdict
self.shape = "nxp,nyp,nzp" # 3차원 그리드 구조
self.dim = ":,:,:"
def depends_on (self, a):
return (a in self.dep)
def is_primary (self):
return self.prime
def not_primary (self):
return not self.prime
def is_fluctuation (self):
return self.fluc
def export_on (self):
return False
def __repr__ (self):
return self.name
def checkFluctuation (self):
"""본 변수 혹은 의존하고 있는 하위 기호들 중에 변동량(Fluctuation) 관련 계산이 개입되어 있는지
상향식으로 전파 추적하는 재귀 메서드입니다.
"""
fset = set([])
for d in map(self.fdict.get, self.dep):
fset.update(d.checkFluctuation())
if self.is_fluctuation() or len(fset) > 0:
fset.add(self.name)
return fset
def depClosure (self):
"""해당 변수를 계산하기 위해 선행 계산되어야 하는 모든 하위 변수 노드들을
재귀적으로 타고 내려가 총합 폐쇄 집합(Closure Set)으로 묶어 반환합니다.
"""
fset = set(self.dep)
for d in self.dep:
fset.update(self.fdict[d].depClosure())
return fset
def code_decl (self):
real_array_decl = "real(real64), allocatable, dimension({1}) :: {0}"
return real_array_decl.format(self.name, self.dim)
def code_alloc (self):
return make_allocate(self.name, self.shape)
def code_free (self):
real_array_free = "deallocate({})"
return real_array_free.format(self.name)
class FieldExporter (object):
"""물리 필드 데이터를 병렬 분산 디스크 시스템으로 직접 추출(Export)하는 고성능 MPI-IO 서브루틴 블록을
생성하는 템플릿 처리 클래스입니다.
1. Subarray 방식 (기본값): MPI_TYPE_CREATE_SUBARRAY를 사용하여 3차원 격자의 전체 물리 도메인에서
현재 MPI 랭크가 맡고 있는 일부분의 3D 블록을 통째로 오프셋과 크기를 지정해 병렬 쓰기함으로써 병목을 최소화합니다.
2. Legacy Copy 방식 (특정 X 좌표들만 샘플링하여 쓸 때): 전체 배열 중 특정 좌표들만 추출하여
임시 전송 배열로 카피한 뒤 단일 스트림 형태로 쓰기를 수행합니다.
"""
mpi_io_decl='''
! field exporter common
integer(kind=MPI_OFFSET_KIND) :: offset
'''
# Subarray 버전 템플릿
fmt_decl_subarray='''
! - file_handles and mpi_infos
integer(kind=MPI_INTEGER_KIND) :: {{ field_name }}_fh
integer(kind=MPI_INTEGER_KIND) :: {{ field_name }}_info
integer(kind=MPI_INTEGER_KIND) :: {{ field_name }}_filetype
'''
fmt_init_subarray='''
! init subarray datatype for {{ field_name }}
block
integer(4) :: sizes(3), subsizes(3), starts(3)
call MPI_INFO_CREATE({{ field_name }}_info, mpi_err)
call MPI_FILE_OPEN(MPI_COMM_TASK,'export-{{ field_name }}.dat',MPI_MODE_WRONLY+MPI_MODE_CREATE,{{ field_name }}_info,{{ field_name }}_fh,mpi_err)
sizes = (/ nxp, nyp, nzp /)
subsizes = (/ {{ len_xpts }}, {{ ye }} - {{ ys }} + 1, {{ ze }} - {{ zs }} + 1 /)
starts = (/ {{ xs }} - 1, {{ ys }} - 1, {{ zs }} - 1 /)
call MPI_TYPE_CREATE_SUBARRAY(3, sizes, subsizes, starts, MPI_ORDER_FORTRAN, MPI_REAL8, {{ field_name }}_filetype, mpi_err)
call MPI_TYPE_COMMIT({{ field_name }}_filetype, mpi_err)
end block
'''
fmt_final_subarray='''
! finalize
call MPI_FILE_CLOSE({{ field_name }}_fh, mpi_err)
call MPI_INFO_FREE({{ field_name }}_info, mpi_err)
call MPI_TYPE_FREE({{ field_name }}_filetype, mpi_err)
'''
fmt_calc_subarray='''
! write to file via MPI Subarray
count = ({{ len_xpts }}) * ({{ ye }} - {{ ys }} + 1) * ({{ ze }} - {{ zs }} + 1)
offset = export_offset(fidx) * count * 8
call MPI_FILE_WRITE_AT({{ field_name }}_fh, offset, {{ work_array }}, 1, {{ field_name }}_filetype, mpi_status, mpi_err)
'''
# Legacy copy 버전 (fallback)
fmt_decl_legacy='''
! - file_handles and mpi_infos
integer(kind=MPI_INTEGER_KIND) :: {{ field_name }}_fh
integer(kind=MPI_INTEGER_KIND) :: {{ field_name }}_info
! - buffer
real(real64), allocatable, dimension(:,:,:) :: {{ field_name }}_export_array
integer, allocatable, dimension(:) :: {{ field_name }}_xpts
'''
fmt_init_legacy='''
! init
call MPI_INFO_CREATE({{ field_name }}_info, mpi_err)
call MPI_FILE_OPEN(MPI_COMM_TASK,'export-{{ field_name }}.dat',MPI_MODE_WRONLY+MPI_MODE_CREATE,{{ field_name }}_info,{{ field_name }}_fh,mpi_err)
allocate({{ field_name }}_export_array(1:{{ len_xpts }},{{ ys }}:{{ ye }},{{ zs }}:{{ ze }}), stat=ierr)
if (ierr /= 0) then
write(0,*) 'Error: allocation of {{ field_name }}_export_array failed on process', myid
call MPI_ABORT(MPI_COMM_TASK, 1, mpi_err)
end if
{{ field_name }}_export_array = 0.
allocate({{ field_name }}_xpts(1:{{ len_xpts }}), stat=ierr)
if (ierr /= 0) then
write(0,*) 'Error: allocation of {{ field_name }}_xpts failed on process', myid
call MPI_ABORT(MPI_COMM_TASK, 1, mpi_err)
end if
{{ xpts_init }}
'''
fmt_final_legacy='''
! finalize
call MPI_FILE_CLOSE({{ field_name }}_fh, mpi_err)
call MPI_INFO_FREE({{ field_name }}_info, mpi_err)
deallocate({{ field_name }}_export_array)
deallocate({{ field_name }}_xpts)
'''
fmt_calc_legacy='''
! copy to array for export
do k = {{ zs }}, {{ ze }}
do j = {{ ys }}, {{ ye }}
do i = 1, {{ len_xpts }}
{{ field_name }}_export_array(i,j,k) = {{ work_array }}({{ field_name }}_xpts(i),j,k)
end do
end do
end do
! write to file
count = ({{ len_xpts }}) * ({{ ye }} - {{ ys }} + 1) * ({{ ze }} - {{ zs }} + 1)
offset = export_offset(fidx) * count * 8
call MPI_FILE_WRITE_AT({{ field_name }}_fh, offset, {{ field_name }}_export_array, count, MPI_REAL8, mpi_status, mpi_err)
'''
def __init__ (self, name, attr, parent):
self.name = name
self.attr = attr
self.parent = parent
self.params = dict(attr)
self.params.setdefault("xs", 1)
self.params.setdefault("xe", "nxp")
self.params.setdefault("ys", 1)
self.params.setdefault("ye", "nyp")
self.params.setdefault("zs", 1)
self.params.setdefault("ze", "nzp")
self.params.setdefault("field_name", self.name)
self.params.setdefault("len_xpts", f"({self.params['xe']} - {self.params['xs']} + 1)")
fmt_xpts_init = '''
do i = {{ xs }}, {{ xe }}
{{ field_name }}_xpts(i-{{ xs }}+1) = i
end do
'''
self.params.setdefault("xpts_init", Template(fmt_xpts_init).render(**self.params))
try:
# Sampling at listed x coordinates
fmt_xpts_init_list = "{{ field_name }}_xpts = (/ {{ list_xpts }} /)"
import sys
raw_xpts = self.params["xpts"]
int_xpts = list(map(int, raw_xpts.split()))
len_xpts = len(int_xpts)
self.params["len_xpts"] = len_xpts
self.params["list_xpts"] = ",".join(map(str, int_xpts))
self.params["xpts_init"] = Template(fmt_xpts_init_list).render(**self.params)
except KeyError:
pass
self.use_subarray = ("xpts" not in self.params)
def code (self):
self.params["work_array"] = self.parent.array
if self.use_subarray:
return Template(FieldExporter.fmt_calc_subarray).render(**self.params)
else:
return Template(FieldExporter.fmt_calc_legacy).render(**self.params)
def code_decl (self):
if self.use_subarray:
return Template(FieldExporter.fmt_decl_subarray).render(**self.params)
else:
return Template(FieldExporter.fmt_decl_legacy).render(**self.params)
def code_alloc (self):
if self.use_subarray:
return Template(FieldExporter.fmt_init_subarray).render(**self.params)
else:
return Template(FieldExporter.fmt_init_legacy).render(**self.params)
def code_free (self):
if self.use_subarray:
return Template(FieldExporter.fmt_final_subarray).render(**self.params)
else:
return Template(FieldExporter.fmt_final_legacy).render(**self.params)
class Field (FieldBase):
"""일반 대입 계산 변수를 관리하며, 3차원 격자점 루프 코드 생성을 담당하는 핵심 클래스입니다.
SymPy 최적화 및 CSE 적용 코드를 루프 본문에 결합합니다.
"""
def __init__ (self, name, attr, exp, fdict):
super(Field,self).__init__(name, fdict)
self.attr = attr
self.exp = exp
self.fluc, self.dep, self.derivs = ExpInspector.inspect(exp)
self.comment = self.name + " = " + ExpToCode(self.fdict).transform(self.exp)
self.latex_given = self.attr.get("latex")
if self.latex_given is None:
self.latex = ExpToLatex(self.fdict).transform(self.exp)
else:
self.latex = self.latex_given
self.exporter = None
try:
if self.attr["export"]:
self.exporter = FieldExporter(self.name, self.attr, self)
except KeyError:
pass
def export_on (self):
return self.exporter is not None
def code (self, alloc=None):
self.array = alloc[self.name] if alloc else self.name
# SymPy 기호 수식 최적화 및 루프 내 CSE 적용 스칼라 추출
opt = SympyOptimizer.get_instance(self.fdict)
rhs, cse_decls, cse_assigns = opt.optimize_field(self.name, alloc)
decls_str = "\n".join(cse_decls) if cse_decls else ""
assigns_str = "\n".join(cse_assigns) if cse_assigns else ""
# 생성될 Fortran 3차원 루프 템플릿
# CSE 최적화에 의해 추출된 로컬 임시 스칼라 할당문(assigns_str)을 3중 루프 i, j, k 본문 내부에서 먼저 연산하도록 삽입하고,
# 이 변수들을 사용하여 우변(rhs)을 연산하고 최종 배열(array(i,j,k))에 저장합니다.
real_array_loop = """
! {{ comment }}
{% if decls_str -%}
block
{{ decls_str | indent(4, True) }}
{%- endif %}
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
{% if assigns_str -%}
{{ assigns_str | indent(4, True) }}
{{ array }}(i,j,k) = {{ rhs }}
{%- else -%}
{{ array }}(i,j,k) = {{ rhs }}
{%- endif %}
end do
end do
end do
{% if decls_str -%}
end block
{%- endif %}
"""
calculation_code = Template(real_array_loop).render(
comment=self.comment,
decls_str=decls_str,
assigns_str=assigns_str,
array=self.array,
rhs=rhs
)
export_code = ( self.exporter.code() if self.export_on() else "")
return calculation_code + export_code
class FluctuationField (FieldBase):
"""물리 필드의 난류 변동 성분(Fluctuation, u' = u - <u_w>)을 계산하기 위한 변수 클래스입니다.
수식 내의 u' 기호를 평균량과의 차이 수식으로 팽창하여 할당합니다.
"""
def __init__ (self, w, field, fset, fdict):
super(FluctuationField,self).__init__(self.id(w,field), fdict)
if w is not None:
self.w = w + "_"
else:
self.w = ""
self.field = fdict[field]
self.dep = self.field.dep - fset
for df in self.field.dep & fset:
self.dep.add(self.id(w,df))
self.comment = ExpToCode(self.fdict).transform(self.field.exp)
if self.field.is_fluctuation():
self.comment = self.comment.format(self.w)
self.comment = self.name + " = " + self.comment
def code (self, alloc=None):
self.array = alloc[self.name] if alloc else self.name
rhs = ExpToCode(self.fdict).transform(self.field.exp)
if self.field.is_fluctuation():
rhs = rhs.format(self.w)
real_array_loop = """
! {{ comment }}
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
{{ array }}(i,j,k) = {{ rhs }}
end do
end do
end do
"""
return Template(real_array_loop).render(
comment=self.comment,
array=self.array,
rhs=rhs
)
@classmethod
def id (cls, w, field):
if w:
name = "{}____{}_avg".format(field, w)
else:
name = "{}____avg".format(field)
return name
class PrimaryField (FieldBase):
"""격자 정보나 외부 물리계 수치 데이터(u, v, w, T 등) 파일에서 사전에 로드하여
메모리에 상주하는 기본 원본 입력 필드 클래스입니다.
자체 계산 루프나 동적 할당 코드를 직접 생성하지 않습니다.
"""
def __init__ (self, name, fdict):
super(PrimaryField,self).__init__(name, fdict)
self.derivs = set([])
self.prime = True
self.latex = name
self.latex_given = None
def code (self, alloc=None):
return "! {} is read from file".format(self.name)
def code_decl (self):
return "! {} is read from file".format(self.name)
def code_alloc (self):
return "! {} is read from file".format(self.name)
def code_free (self):
return "! {} is read from file".format(self.name)
class DerivedField (FieldBase):
"""수치 공간 미분(ddx, d2dy 등)을 수행하여 계산되는 유도 필드 클래스입니다.
Fortran 수치 차분 패키지 서브루틴(Compact.f90 에 구현된 dfnonp, dfp 등)의
동적 호출 코드를 출력합니다.
"""
def __init__ (self, op, v, fdict):
name = "{}_{}".format(op, v)
super(DerivedField,self).__init__(name, fdict)
self.op = op
self.v = v
self.dep = set([v])
fmt = r"\partial_{{{}}}"
coord = op[-1]
partial = fmt.format(coord + coord if len(op) > 3 else coord)
self.latex = partial + "(" + fdict[v].latex + ")"
def code (self, alloc=None):
self.array = alloc[self.name] if alloc else self.name
varray = alloc[self.v] if alloc else self.v
# 예: call ddx ( xyzbuffer0, u ) 형태의 서브루틴 호출 코드를 작성
return "call {0} ( {2}, {1} )".format(self.op, varray, self.array)
class AveragedField (FieldBase):
"""특정 물리 필드를 격자의 동질 차원(예: X 방향 1D선상 평균)에 대해
공간 통계 평균(Average) 연산을 수행하는 1차원 필드 클래스입니다.
"""
@classmethod
def id (cls, w, tgt):
if w:
return "{}_avg_{}".format(w, tgt)
else:
return "avg_{}".format(tgt)
def __init__ (self, w, tgt, fdict):
name = self.id(w,tgt)
super(AveragedField,self).__init__(name, fdict)
self.shape = "nxp" # Y, Z 차원을 평균화하여 날려버리므로 X방향 크기인 nxp 1차원 배열이 됨
self.dim = ":"
self.target = tgt
tfield = fdict[tgt]
self.fset = tfield.checkFluctuation()
self.latex = r"\left\langle {} \right\rangle".format(tfield.latex)
if not self.fset:
self.tgt = tgt
self.dep.add(tgt)
else:
ftgt = FluctuationField.id(w,tgt)
self.tgt = ftgt
self.dep.add(ftgt)
self.weighted = w
if w:
self.w = fdict[w]
self.dep.add(w)
self.latex += ("_{{{}}}".format(w))
def code (self, alloc=None):
# 3차원 루프를 돌며 Y, Z 축에 대해 값을 누적합 연산하는 구문
avg_array_sum = """
do k = 1, nzp
do j = 1, nyp
do i = 1, nxp
{{ name }}(i) = {{ name }}(i) + {{ arrname }}
end do
end do
end do
"""
arrname = self.fdict[self.tgt].array + "(i,j,k)"
if self.weighted is not None:
arrname = arrname + " * " + self.w.array + "(i,j,k)"
return Template(avg_array_sum).render(name=self.name, arrname=arrname)
def code_avg (self):
# 병렬 랭크 간 MPI_ALLREDUCE를 통해 X방향 라인별 총합을 싱크한 뒤,
# 전체 그리드 면적(denum = nyp * nzp) 및 가중 평균 값으로 나누어 실제 수학적 평균값 산출
avg_array_divide = """
call MPI_ALLREDUCE(MPI_IN_PLACE, {{ name }}, nxp, MPI_REAL8, MPI_SUM, MPI_COMM_TASK, mpi_err)
{{ name }} = {{ name }} {{ dWeight }} / denum
"""
dWeight = (f"/ avg_{self.weighted}" if self.weighted else "")
return Template(avg_array_divide).render(name=self.name, dWeight=dWeight)
def isWeighted (self):
return self.weighted is not None
def pass1 (self):
return not self.pass2()
def pass2 (self):
# 수식 내에 u' 등 변동량 항이 포함되어 있으면,
# 원본 u의 평균치(<u_w>) 계산이 완료된 '이후'에만 연산이 가능하므로 Pass 2 대상이 됨
return len(self.fset) > 0
class Stage1():
"""컴파일러 1단계: DSL 파서 AST 트리로부터 필드 변수 정의 정보를 추출 및 적재합니다.
Lark 구문 파서가 생성한 추상 구문 트리(AST)를 CollectDefinitions Visitor 클래스로 횡단 스캔하여,
기본 변수(Primary), 수식 유도 변수(Derived), 평균 대상 사양(Averaged)을 모아 구조체 딕셔너리로 저장합니다.
"""
def __init__ (self, raw_tree):
"""Stage 1 초기화 및 AST 정의 스캔을 시작합니다.
Args:
raw_tree (lark.Tree): Lark 구문 파서가 해석 완료한 DSL 추상 구문 트리.
"""
self.primary = set([])
self.derived = {}
self.averaged = {}
# CollectDefinitions Visitor를 구동하여 AST 트리 내의 정의 스캔 시작
CollectDefinitions( self.primary, self.derived, self.averaged ).visit(raw_tree)
def __repr__ (self):
return "\n".join(map(str, [self.primary, self.derived, self.averaged]))
class Stage2():
"""컴파일러 2단계: 수식 내의 미분 및 변동 연산자를 감지하여 해당 물리 필드로 확장 처리합니다.
1단계에서 추출된 변수 식 내부에 ddx(u) 와 같은 수치 미분항이 검출되면, 이를 연산할 중간 미분 필드(DerivedField)를
자동 생성하고 의존성 체인에 등록합니다. 또한, u'와 같은 변동 성분이 포함되어 있으면
물리량에서 평균치를 뺀 변동량 계산용 임시 필드(FluctuationField)로 확장 생성하여 등록합니다.
"""
def __init__ (self, src):
"""Stage 2 초기화 및 변수 수식 확장을 가동합니다.
Args:
src (Stage1): 완료된 Stage 1의 AST 데이터 세트.
"""
self.src = src
self.primary = src.primary
self.derived = src.derived
self.derivative = {}
self.averaged = {}
# 1. 계산식 내에 존재하는 고차 차분 미분항(ddx, d2dy 등)을 찾아 중간 DerivedField로 등록
dset = set([])
for k, v in self.derived.items():
dset.update(v.derivs)
for tup in dset:
a = DerivedField(tup[0], tup[1], self.derived)
self.derived[a.name] = a
self.derivative[tup] = a
# 2. 통계 물리량 평균화 대상 변수들을 AveragedField 구조체로 구성하고,
# 변동량 계산이 필요한 경우 FluctuationField로 등록
for w, tgts in src.averaged.items():
for t in tgts:
a = AveragedField(w, t, self.derived)
self.averaged[a.name] = a
# 평균 편차가 동반된 항들에 대해 FluctuationField 생성
for ff in a.fset:
b = FluctuationField(w, ff, a.fset, self.derived)
self.derived[b.name] = b
def __repr__ (self):
return "\n".join(map(str, [self.derived, self.derivative, self.averaged]))
def dependency (self):
"""프로그램 내 모든 필드 간의 1차 의존 관계 그래프(Dependency Graph)를 추출하여 딕셔너리로 반환합니다."""
dgraph = {}
for k,v in self.derived.items():
dgraph[k] = v.dep
for k,v in self.averaged.items():
dgraph[k] = v.dep
return dgraph
class Stage3():
"""컴파일러 3단계: 연산 간의 의존성을 해결하여 3D 계산 루프의 최적 실행 순서를 계산합니다.
1. 루프 분할(Loop Splitting): 난류 평균량 연산의 인과관계에 따라 물리 연산을 두 개의 단계(Pass)로 쪼갭니다.
- Pass 1 (평균 계산 전): 입력 물리량들로부터 1차 수식을 연산하고, 1차 평균값(예: <u>)을 누적합하는 연산들.
- Pass 2 (평균 계산 후): 구해진 평균값을 활용하여 변동량(u' = u - <u>)을 계산하고,
이를 통반한 최종 결합 물리량을 추가적으로 계산/평균화하는 연산들.
2. 최적화 위상 정렬(Topological Sort): 각 Pass 내에서 변수의 선행 의존 관계가 완전히 해결된
올바른 순서로 순차 연산 코드가 생성되도록 정렬 알고리즘을 가동합니다.
"""
def __init__ (self, src):
"""Stage 3 초기화 및 위상 정렬을 수행하여 Pass 1, 2 실행 순서를 결정합니다.
Args:
src (Stage2): 확장된 변수 테이블을 가진 Stage 2 결과 객체.
"""
self.src = src
self.primary = src.primary
self.derived = src.derived
self.averaged = src.averaged
self.dependency = src.dependency()
# 이름 충돌 방지 검증
assert set(self.derived.keys()).isdisjoint(self.averaged.keys())
# Pass 1과 Pass 2 대상 평균화 변수 논리 분할
self.avg1 = set(filter(AveragedField.pass1, self.averaged.values()))
self.avg2 = set(filter(AveragedField.pass2, self.averaged.values()))
# Pass 1 위상 정렬: Pass 1 평균 변수들의 연산에 관여하는 모든 종속 관계를 수집하여 정렬
pass1calc = set(map(repr, self.avg1))
for x in self.avg1:
pass1calc.update(x.depClosure())
self.pass1 = self.sort_vars_new(self.dependency, pass1calc - self.primary)
# Pass 2 위상 정렬: Pass 2 평균 변수(변동 연산 연계)들의 연산에 관여하는 종속성을 정렬
pass2calc = set(map(repr, self.avg2))
for x in self.avg2:
pass2calc.update(x.depClosure())
self.pass2 = self.sort_vars_new(self.dependency, pass2calc - self.primary)
def __repr__ (self):
return "\n".join(map(str, [self.pass1, self.pass2]))
def sort_vars (self, dependency, group):
"""일반적인 단순 위상 정렬 알고리즘입니다."""
order = []
remain = list(group)
remain.sort()
while len(remain) > 0:
for v in remain:
if dependency[v].isdisjoint(remain):
order.append(v)
remain.remove(v)
return order
def calc_size (self, ordered, remaining):
"""위상 정렬 도중 활성화되는 임시 라이브 변수들의 상대적 가중치 스케일을 계산합니다.
메모리 상주 대역을 최소화하는 최적의 정렬 후보를 결정할 때 척도로 이용됩니다.
"""
count = 0
dep_union = set([])
g = self.dependency
for v in remaining:
dep_union |= set(g[v])
for v in ordered:
if v in dep_union:
count += 1
return count
def sort_vars_new (self, dependency, group):
"""메모리에 동시에 상주하여 참조되어야 하는 활성 라이브 변수의 개수를
최소화하도록 특수 설계된 메모리 인지형 위상 정렬(Memory-Aware Topological Sort) 알고리즘입니다.
매 정렬 단계마다 현재까지 결정된 실행 흐름(ordered)과 향후 연산되어야 할 대기 변수(remaining) 사이의
의존도 영향 범위(calc_size)를 계산하여, 생명 주기 창(Liveness Window)이 짧아지도록 후보군 중 최선의 변수를
선택해 나가며 실행 순서 목록을 빌드합니다.
"""
order = []
remain = list(group)
remain.sort()
while len(remain) > 0:
candidate = []
for v in remain:
if set(dependency[v]).isdisjoint(remain):
candidate.append(v)
impact = {}
size0 = self.calc_size(set(order), set(remain))
for v in candidate:
impact[v] = self.calc_size(set(order) | set([v]), set(remain) - set([v])) - size0
candidate.sort(key=impact.get)
# 활성 변수 상주 가중치 영향도가 가장 적은 후보를 우선 배치
order.append(candidate[0])
remain.remove(candidate[0])
return order
def print_program (self):
"""Stage 3 결과를 바탕으로, 최적화가 적용되지 않은 평이한 단일 3D 배열식으로
모든 중간 임시 변수를 개별 동적 할당하여 계산하는 원시 Fortran 소스코드를 인쇄합니다.
"""
allvar = dict(self.derived)
allvar.update(self.averaged)
decl = "\n".join(allvar[v].code_decl() for v in set(self.pass1+self.pass2))
alloc = "\n".join(allvar[v].code_alloc() for v in set(self.pass1+self.pass2))
free = "\n".join(allvar[v].code_free() for v in set(self.pass1+self.pass2))
calc1 = "\n".join(allvar[v].code() for v in self.pass1)
calc2 = "\n".join(allvar[v].code() for v in self.pass2)
set1 = [a.name for a in filter(AveragedField.pass1, self.averaged.values())]
set2 = [a.name for a in filter(AveragedField.pass2, self.averaged.values())]
set1.sort()
set2.sort()
avg1 = "\n".join(allvar[v].code_avg() for v in set1)
avg2 = "\n".join(allvar[v].code_avg() for v in set2)
hfmt = 'character (len = *), parameter :: output_header="{}"'
declh = hfmt.format(" ".join(["x"] + set1 + set2))
avg_array_write = '''
integer :: i
open (200, file="qEdge_X.dat")
write (200,*) output_header
do i=1,nxp
write (200,'({0}e20.10)') real(i)*hxp, {1}
end do
close (200)
'''
avgarr = "{}(i)"
write_avg = avg_array_write.format(
len(self.averaged)+1,
", ".join(map(avgarr.format, set1+set2))
)
md = {}
md["module_name"] = "terms"
md["module_data"] = "\n".join((declh, decl))
md["module_init"] = alloc
md["module_finalize"] = free
md["module_pass1"] = calc1
md["module_pass1_avg"] = avg1
md["module_pass2"] = calc2
md["module_pass2_avg"] = avg2
md["module_write_result"] = write_avg
return md
class Stage4():
"""컴파일러 4단계 (최종): 메모리 Pooling 버퍼 할당 및 SymPy 최적화 코드를 최종 결합하여 컴파일을 완료합니다.
[핵심 최적화 기술]
1. SymPy 수식 정리 및 의존 기호 교정:
- SymPy 기호 엔진을 통해 최적화 대상 필드들의 수식 트리를 전면 간소화합니다.
- 기호화 과정에서 대입 치환(substitution)으로 바뀐 수식들의 의존성을 재추적하여
컴파일러의 dependency 그래프를 정확한 상태로 최신화합니다.
2. 생명 주기 분석 (Liveness Analysis):
- 3D 배열은 크기가 막대하여 루프 중간에 생성되는 수십 개의 임시 변수들에게 각자 별도 배열을 할당하면
프로그램 실행이 불가능해집니다.
- 임시 변수가 할당되어 계산되는 시작점과 그 변수를 다른 수식들이 완전히 가져다 쓴 뒤
종료되는 참조 종점(Liveness window)을 2차원 부울 마스크 이미지 형태로 추적합니다.
3. 메모리 Pooling 및 버퍼 할당 (Array Register Allocation):
- Liveness 마스크 분석을 통해 동일 시점에 공존하여 메모리에 살아있어야 하는 임시 변수의 최대 중첩 개수를 파악하여
필요한 버퍼의 총수(self.narr)를 계산합니다.
- 재사용 가능한 버퍼 명칭 풀(pool)을 정의하고, Liveness에 의해 변수의 생명이 다하면 사용하던 버퍼명을 풀에 반환하며,
새로운 변수가 활성화될 때 풀에서 빈 버퍼명을 분배(Pop)하는 방식으로 3차원 물리 필드 메모리 오버헤드를
O(1) 수준의 버퍼 몇 개 세트로 극도로 낮춥니다 (xyzbuffer0, xyzbuffer1 등).
"""
def __init__ (self, src):
self.src = src
self.primary = src.primary
self.derived = src.derived
self.averaged = src.averaged
self.dependency = src.dependency
self.avg1 = src.avg1
self.avg2 = src.avg2
self.pass1 = src.pass1
self.pass2 = src.pass2
# 1. SymPy 수식 최적화 엔진 초기화
opt = SympyOptimizer.get_instance(self.derived)
opt.set_averaged(self.averaged)
# 2. SymPy가 대입식 치환(Substitution) 과정에서 제거한 불필요한 의존성 관계를
# 의존성 그래프에 즉각 반영하여 실제 계산을 위한 의존성 체인을 슬림하게 정리합니다.
updated_dependency = {}
for name, dep_set in self.dependency.items():
if name in self.derived and isinstance(self.derived[name], Field):
expr = opt.get_sympy_expr(name)
# 실제 정리된 SymPy 식에 잔존하는 자유 기호 명칭들만 추출
free_sym_names = {sym.name for sym in expr.free_symbols}
valid_deps = {dep for dep in free_sym_names if dep in self.derived or dep in self.primary}
updated_dependency[name] = valid_deps
else:
updated_dependency[name] = dep_set
self.dependency = updated_dependency
self.array_name = "xyzbuffer{}" # Pooling용 공용 3차원 버퍼 이름 포맷
# 3. Pass 1 및 Pass 2 연산 순서 배열들에 대해 각각 버퍼 공유 매핑(Pooling) 수행
narr1, alloc1 = (self.allocate_arr(self.pass1))
narr2, alloc2 = (self.allocate_arr(self.pass2))
# 전체 프로그램에서 필요한 동적 공유 3D 버퍼 배열의 최대 크기 설정
self.narr = max(narr1, narr2)
self.alloc1 = alloc1 # Pass 1 변수명 -> 공용 xyzbuffer 명칭 매핑 정보
self.alloc2 = alloc2 # Pass 2 변수명 -> 공용 xyzbuffer 명칭 매핑 정보
def liveness (self, l1, g):
"""위상 정렬된 변수 시퀀스 `l1`과 의존성 그래프 `g`를 이용하여
각 변수가 루프 타임라인 내내 살아있어야 하는 구간(Liveness window)을 2D 부울 행렬로 구성합니다.
행렬의 (i, j)가 True이면, i번째 변수가 j번째 연산 시점에 살아있어야 함을 의미합니다.
"""
import numpy as np
img = np.zeros((len(l1), len(l1)))
for i, v in enumerate(l1):
for j in range(i, len(l1)):
# 만약 i번째 변수가 훗날 j번째 변수를 계산할 때의 하위 의존성(g[l1[j]])에 포함되어 있으면,
# i번째 변수는 j번째 연산 시점까지 생명이 연장되어 살아있어야(live) 합니다.
img[i,i:j+1] = img[i,i:j+1] + (1 if v in g[l1[j]] else 0)
return img > 0
def allocate_arr (self, l):
"""Liveness 분석 데이터를 바탕으로 동적 3D 배열 레지스터 할당(Pooling)을 수행합니다.
동시 활성 상태를 파악하여 최소 버퍼 개수(narr)를 도출하고, 변수 간 버퍼 공유 관계(var2arr)를 리턴합니다.
"""
import numpy as np
dg = self.dependency
# 1. 생명 주기 마스크 행렬 획득
mask = self.liveness(l, dg)
# 2. 열 방향 합계를 내어 동시에 생존해야 하는 임시 변수의 최대 개수(동시 상주 배열 피크치)를 산출
try:
narr = mask.astype(int).sum(axis=0).max()
except ValueError:
narr = 0
# 사용 가능한 공유 버퍼 풀 정의
array_pool = set([self.array_name.format(i) for i in range(narr)])
# 각 계산 단계별로 살아있는 변수들의 집합 구성
livesets = [set([])] + [set(np.asarray(l)[row]) for row in mask.T]
# 기본 입력 변수들은 파일로부터 직접 매핑되므로 원래 명칭 유지
var2arr = { p : p for p in self.primary }
# 루프 연산 순서 타임라인을 한 단계씩 밟아가며 동적 버퍼 매핑 시뮬레이션 수행
for i, (s0, s1) in enumerate(zip(livesets[:-1], livesets[1:])):
# 이전 단계(s0)에는 살아있었으나 이번 단계(s1)에서는 죽은 변수들(s0 - s1)의
# 할당 버퍼를 다시 사용 가능한 array_pool로 회수(반환)합니다.
array_pool.update(map(var2arr.get, s0 - s1))
# 이번 단계에 새로 생성되는 변수들(s1 - s0)에게 풀에서 노는 버퍼명을 Pop하여 할당(배분)합니다.
for new in s1 - s0:
var2arr[new] = array_pool.pop()
return narr, var2arr
def work_array_codes (self):
"""생성된 Fortran 모듈 상단에 선언 및 할당/해제될 공용 3D 버퍼 배열(xyzbuffer0, ...)의
선언문, allocate문, deallocate문 코드를 작성합니다.
"""
array_names = [self.array_name.format(i) for i in range(self.narr)]
real_array_decl = "real(real64), allocatable, dimension(:,:,:) :: {0}"
decl = "\n".join([real_array_decl.format(v) for v in array_names])
alloc = "\n".join([make_allocate(v, "nxp,nyp,nzp") for v in array_names])
free = "\n".join(["deallocate({})".format(v) for v in array_names])
return decl, alloc, free
def write_avg_codes (self, avglist):
"""연산 최종 완료 후 평균화가 끝난 1D 결과 데이터셋 및 미분량을
텍스트 파일(qEdge_X.dat, d1.dat, d2.dat)에 형식에 맞추어 라이팅하는 Fortran 서브루틴 블록을 생성합니다.
"""
avg_array_write = '''
real(real64), dimension(nxp) :: xbuffer
integer :: i
open (200, file="qEdge_X.dat")
write (200,*) output_header
do i=1,nxp
write (200,'({{ num_args }}e20.10)') real(i)*hxp, {{ formatted_avglist }}
end do
close (200)
open (200, file="d1.dat")
{{ deriv1_lines }}
close (200)
open (200, file="d2.dat")
{{ deriv2_lines }}
close (200)
'''
avgarr = "{}(i)"
deriv1_avgarr = """call ddx1d ( xbuffer, {} ) ; write (200,*) xbuffer"""
deriv2_avgarr = """call d2dx1d ( xbuffer, {} ) ; write (200,*) xbuffer"""
num_args = len(self.averaged) + 1
formatted_avglist = ", ".join(map(avgarr.format, avglist))
deriv1_lines = "\n".join(map(deriv1_avgarr.format, avglist))
deriv2_lines = "\n".join(map(deriv2_avgarr.format, avglist))
write_avg = Template(avg_array_write).render(
num_args=num_args,
formatted_avglist=formatted_avglist,
deriv1_lines=deriv1_lines,
deriv2_lines=deriv2_lines
)
return write_avg
def print_program (self):
"""4단계 최적화(공유 Pooling 버퍼 매핑 및 SymPy CSE 치환)가 완료된
가장 고성능의 최종 Fortran 모듈 코드를 생성하기 위해 Jinja2용 딕셔너리를 빌드합니다.
"""
opt = SympyOptimizer.get_instance(self.derived)
opt.set_averaged(self.averaged)
allvar = dict(self.derived)
allvar.update(self.averaged)
# 평균 변수들을 순서대로 분배
set1 = sorted([a.name for a in filter(AveragedField.pass1, self.averaged.values())])
set2 = sorted([a.name for a in filter(AveragedField.pass2, self.averaged.values())])
# 외부 디스크 파일 익스포트 활성화 여부 확인
set_export_on = list(filter(lambda x: x.export_on(), self.derived.values()))
ffmt = 'logical, parameter :: pass2_required={}'
declf = ffmt.format('.true.' if len(set2) > 0 else '.false.')
hfmt = 'character (len = *), parameter :: output_header="{}"'
declh = hfmt.format(" ".join(["x"] + set1 + set2))
# 공용 Pooling 3D 버퍼 배열의 선언/할당 코드 생성
declarr, allocarr, freearr = self.work_array_codes()
# 1차원 평균 물리량 배열들의 선언/할당 코드 생성
declavg = "\n".join(self.averaged[v].code_decl() for v in sorted(self.averaged))
allocavg = "\n".join(self.averaged[v].code_alloc() for v in sorted(self.averaged))
freeavg = "\n".join(self.averaged[v].code_free() for v in sorted(self.averaged))
# 병렬 파일 쓰기(MPI Subarray)를 위한 MPI 리소스 선언/할당 코드 생성
decl_export = "\n".join(v.exporter.code_decl() for v in set_export_on)
alloc_export = "\n".join(v.exporter.code_alloc() for v in set_export_on)
free_export = "\n".join(v.exporter.code_free() for v in set_export_on)
# Pass 1과 Pass 2 루프 내부 본문 연산 코드들을 생성 (각자 버퍼 맵 alloc1, alloc2 적용)
sub_calc1 = "\n".join(allvar[v].code(self.alloc1) for v in self.pass1 if v in self.averaged or v in self.alloc1)
sub_calc2 = "\n".join(allvar[v].code(self.alloc2) for v in self.pass2 if v in self.averaged or v in self.alloc2)
# 평균 누적 연산 코드 생성
sub_avg1 = "\n".join(allvar[v].code_avg() for v in set1)
sub_avg2 = "\n".join(allvar[v].code_avg() for v in set2)
# 파일 쓰기 루틴 코드 생성
sub_write_avg = self.write_avg_codes(set1+set2)
md = {}
md["module_name"] = "terms"
md["module_data"] = "\n".join((declf, declh, declavg, FieldExporter.mpi_io_decl, decl_export, declarr))
md["module_init"] = "\n".join((allocavg, alloc_export, allocarr))
md["module_finalize"] = "\n".join((freeavg, free_export, freearr))
md["module_pass1"] = sub_calc1
md["module_pass1_avg"] = sub_avg1
md["module_pass2"] = sub_calc2
md["module_pass2_avg"] = sub_avg2
md["module_write_result"] = sub_write_avg
return md
def save_ir (self):
"""코드 생성 디버깅을 위해 중간 컴파일 단계 변수 목록 및 의존성 맵을 JSON 형태 파일로 덤프합니다."""
import json
dg = {k:list(v) for k,v in self.dependency.items()}
with open("ir2.py", "w") as irf:
print("g = ", json.dumps(dg, indent=4), file=irf)
print("l1 = ", json.dumps(self.pass1, indent=4), file=irf)
print("l2 = ", json.dumps(self.pass2, indent=4), file=irf)
print("avg1 = ", json.dumps(list(map(repr,self.avg1)), indent=4), file=irf)
print("avg2 = ", json.dumps(list(map(repr,self.avg2)), indent=4), file=irf)
class Stage5():
"""예비 단계: 확장된 루프 분할 및 연산 재정렬을 수행하는 플레이스홀더 클래스입니다."""
pass
class Stage6():
"""예비 단계: 최종 다이렉트 Fortran 기계 최적화 코드 출력 플레이스홀더 클래스입니다."""
pass