""" 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)하도록 할당하여 메모리 사용량을 최소화합니다. """ import sys 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) 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_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' 기호를 평균량과의 차이 수식으로 팽창하여 할당합니다. """ 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의 평균치() 계산이 완료된 '이후'에만 연산이 가능하므로 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차 평균값(예: ``)을 누적합하는 연산들. - Pass 2 (평균 계산 후): 구해진 평균값을 활용하여 변동량(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)