feat(cryotank): 水素タンク FEM 生成の着手 (DOE + Abaqus 雛形)

docs/HYDROGEN_TANK_SHM.md

@@ -119,11 +119,17 @@ GW グラフ（センサ=ノード）は `build_gw_graph.py` のスキーマを
 ## 7. 1サンプル FEM 実行計画 (One-Sample Plan)
 
 プロジェクト流儀に従い、**まず1サンプルで検証**してからバッチ化する。
+着手済みコード: `src/generate_cryotank_doe.py`（DOE, Abaqus 非依存・テスト済み）、
+`src/generate_cryotank_dataset.py`（Abaqus CAE 雛形, `--dry-run` は非 Abaqus で解析プラン確認可）。
 
 ```bash
-# 1) 生成スクリプト雛形（新規, generate_fairing_dataset.py ベース）
-#    → src/generate_cryotank_dataset.py（healthy 1件）
-abaqus cae noGUI=src/generate_cryotank_dataset.py -- --config doe_cryotank_1sample.json
+# 0) DOE 生成（healthy 1件, Abaqus 不要）
+python src/generate_cryotank_doe.py --healthy_only --n_samples 1 --output doe_cryotank_1sample.json
+# ローカルで解析プランを確認（Abaqus 不要）
+python src/generate_cryotank_dataset.py --dry-run --defect doe_cryotank_1sample.json
+
+# 1) FEM 生成（クラスタ, Abaqus）— 雛形の M1/M2 TODO を実装後に実行
+abaqus cae noGUI=src/generate_cryotank_dataset.py -- --job CryoTank_Healthy_0000 --defect doe_cryotank_1sample.json
 
 # 2) ODB 抽出（既存を流用）
 abaqus python src/extract_odb_results.py --odb abaqus_work/Job-CryoTank-Healthy.odb
@@ -157,8 +163,9 @@ cd src && python build_graph.py --data_dir ../dataset_cryotank_1sample
 
 ## 9. 次アクション
 
+- [x] `src/generate_cryotank_doe.py` / `src/generate_cryotank_dataset.py` 雛形作成（DOE + dry-run はテスト済み）
 - [ ] 本メモのレビュー・数値の確定（材料・寸法・荷重）
-- [ ] `src/generate_cryotank_dataset.py` 雛形作成（healthy 1件）
+- [ ] Abaqus 雛形の M1/M2 TODO 実装（sector part・材料・内圧+熱・欠陥）→ クラスタで healthy 1件生成
 - [ ] M1 検証 → 本メモに結果追記
 - [ ] `ROADMAP.md` / `CLAUDE.md` の Research Lines に「Hydrogen Tank SHM」を正式追加
 

src/generate_cryotank_dataset.py

@@ -0,0 +1,154 @@
+# -*- coding: utf-8 -*-
+# generate_cryotank_dataset.py
+# Abaqus CAE scaffold — H3 cryogenic LH2 tank (Al-Li) FEM with SHM defects.
+#
+# Design note: docs/HYDROGEN_TANK_SHM.md
+# This is the SCAFFOLD for milestone M1 (memo §7): generate ONE healthy sample,
+# eyeball the graph, then batch. Abaqus model-building bodies are marked TODO and
+# should reuse the fairing generator logic (generate_fairing_dataset.py:
+# sector part, symmetry BCs, mesh control; generate_cohesive_fairing.py /
+# generate_czm_sector12.py: cohesive/CZM interfaces for debond & weld flaw).
+#
+# The Abaqus import is guarded so this module can be imported (and --dry-run
+# executed) WITHOUT Abaqus — the cluster runs the real build:
+#   abaqus cae noGUI=src/generate_cryotank_dataset.py -- --job <name> --defect <json>
+# Local plan check (no Abaqus):
+#   python src/generate_cryotank_dataset.py --dry-run --defect doe_cryotank_1sample.json
+
+import argparse
+import json
+import sys
+
+try:
+    from abaqus import *              # noqa: F401,F403
+    from abaqusConstants import *     # noqa: F401,F403
+    from caeModules import *          # noqa: F401,F403
+    from driverUtils import executeOnCaeStartup
+    _HAVE_ABAQUS = True
+except Exception:
+    _HAVE_ABAQUS = False
+
+# ==============================================================================
+# PARAMETERS — Al-Li LH2 tank (design assumption; memo §1-§2, dimensions non-public)
+# ==============================================================================
+RADIUS = 2600.0        # mm, φ ≈ 5.2 m barrel
+H_BARREL = 4000.0      # mm, barrel section modeled
+SECTOR_DEG = 30.0      # deg, symmetric sector (1/12)
+T_WALL = 4.0           # mm, wall thickness (design assumption)
+
+# Al-Li material — RT baseline + cryogenic shift (memo §2; DESIGN ASSUMPTION, verify).
+MAT_ALLI = {
+    "name": "AlLi_2219",
+    "rt": {"E_MPa": 75000.0, "nu": 0.33, "yield_MPa": 400.0, "alpha_per_C": 23.0e-6,
+           "density_t_mm3": 2.70e-9},
+    # cryogenic (20 K) multipliers relative to RT
+    "cryo_factor": {"E": 1.10, "yield": 1.20, "alpha": 0.55},
+}
+
+REFERENCE_TEMP_K = 293.0
+
+
+def material_at(temperature_K):
+    """Return effective Al-Li properties at a temperature (RT interpolation endpoints)."""
+    rt = MAT_ALLI["rt"]
+    if temperature_K <= 20.0:
+        cf = MAT_ALLI["cryo_factor"]
+        return {"E_MPa": rt["E_MPa"] * cf["E"], "nu": rt["nu"],
+                "yield_MPa": rt["yield_MPa"] * cf["yield"],
+                "alpha_per_C": rt["alpha_per_C"] * cf["alpha"],
+                "density_t_mm3": rt["density_t_mm3"], "temperature_K": temperature_K}
+    return {"E_MPa": rt["E_MPa"], "nu": rt["nu"], "yield_MPa": rt["yield_MPa"],
+            "alpha_per_C": rt["alpha_per_C"], "density_t_mm3": rt["density_t_mm3"],
+            "temperature_K": temperature_K}
+
+
+def resolve_plan(defect_params):
+    """Build a human-readable analysis plan from defect_params (no Abaqus needed)."""
+    dtype = defect_params.get("defect_type", "healthy")
+    temp = float(defect_params.get("temperature_K", REFERENCE_TEMP_K))
+    press = float(defect_params.get("pressure_MPa", 0.3))
+    mat = material_at(temp)
+    plan = {
+        "geometry": {"radius_mm": RADIUS, "h_barrel_mm": H_BARREL,
+                     "sector_deg": SECTOR_DEG, "wall_mm": T_WALL},
+        "material": mat,
+        "loading": {"internal_pressure_MPa": press,
+                    "thermal": {"reference_K": REFERENCE_TEMP_K, "operating_K": temp,
+                                "delta_K": temp - REFERENCE_TEMP_K}},
+        "defect": {"type": dtype},
+    }
+    if dtype != "healthy":
+        plan["defect"].update({
+            "theta_deg": defect_params.get("theta_deg"),
+            "z_center_mm": defect_params.get("z_center"),
+            "radius_mm": defect_params.get("radius"),
+            "model": {
+                "weld_flaw": "cohesive/element weakening on weld line",
+                "microcrack": "local stiffness loss + partial contact discontinuity",
+                "insulation_debond": "interface cohesive degradation (reuse CZM)",
+            }.get(dtype, "TODO"),
+        })
+    return plan
+
+
+# ==============================================================================
+# Abaqus build (cluster only) — bodies TODO, reuse fairing generator logic
+# ==============================================================================
+def build_model(job_name, defect_params):
+    if not _HAVE_ABAQUS:
+        raise RuntimeError("Abaqus not available. Use --dry-run for a local plan check, "
+                           "or run via: abaqus cae noGUI=generate_cryotank_dataset.py -- ...")
+    executeOnCaeStartup()
+    plan = resolve_plan(defect_params)
+    # TODO(M1): create sector shell part (reuse generate_fairing_dataset.py sector logic)
+    # TODO(M1): assign Al-Li material from plan["material"]; set expansion (alpha) + ref temp
+    # TODO(M1): sections + assembly + symmetry BCs on circumferential edges
+    # TODO(M1): Step-1 static: internal pressure (plan["loading"]) + Predefined Field
+    #           temperature (operating_K) for cryogenic thermal stress
+    # TODO(M2): insert defect per plan["defect"]["model"]
+    #           - insulation_debond / weld_flaw: cohesive interface (generate_cohesive_fairing.py,
+    #             generate_czm_sector12.py)
+    #           - microcrack: local stiffness reduction + seam
+    # TODO(M1): mesh (continuum/solid shell) + element type
+    # TODO(M1): create + write Job(job_name)
+    raise NotImplementedError("Abaqus model build is a scaffold — implement M1 bodies.")
+
+
+def main():
+    p = argparse.ArgumentParser(description="Cryogenic tank FEM generator (Abaqus)")
+    p.add_argument("--job", type=str, default="CryoTank_Healthy_0000")
+    p.add_argument("--defect", type=str, default=None,
+                   help="JSON: single defect_params, or a DOE {'samples':[...]} (uses sample 0)")
+    p.add_argument("--dry-run", action="store_true", default=False,
+                   help="Print the resolved analysis plan without Abaqus")
+    # Abaqus passes script args after '--'; argparse needs them isolated.
+    argv = sys.argv[1:]
+    if "--" in argv:
+        argv = argv[argv.index("--") + 1:]
+    args = p.parse_args(argv)
+
+    defect_params = {"defect_type": "healthy", "temperature_K": 20.0, "pressure_MPa": 0.3}
+    if args.defect:
+        with open(args.defect) as f:
+            data = json.load(f)
+        if isinstance(data, dict) and "samples" in data:
+            defect_params = data["samples"][0]["defect_params"]
+        elif isinstance(data, dict) and "defect_params" in data:
+            defect_params = data["defect_params"]
+        else:
+            defect_params = data
+
+    if args.dry_run or not _HAVE_ABAQUS:
+        plan = resolve_plan(defect_params)
+        print(json.dumps({"job": args.job, "abaqus_available": _HAVE_ABAQUS, "plan": plan},
+                         indent=2))
+        if not args.dry_run and not _HAVE_ABAQUS:
+            print("\n[note] Abaqus not found — printed plan only. "
+                  "Run on the cluster to build the model.", file=sys.stderr)
+        return
+
+    build_model(args.job, defect_params)
+
+
+if __name__ == "__main__":
+    main()

src/generate_cryotank_doe.py

@@ -0,0 +1,155 @@
+#!/usr/bin/env python3
+# -*- coding: utf-8 -*-
+"""
+Cryogenic Hydrogen-Tank DOE Generation — LHS over defect + operating conditions.
+
+Design note: docs/HYDROGEN_TANK_SHM.md (geometry, materials, loading, defects).
+This is the *first* code step of the hydrogen-tank SHM line (memo §7, milestone M1).
+
+Pure-Python (no Abaqus) — runnable/verifiable outside the cluster. Produces a JSON
+DOE consumed by src/generate_cryotank_dataset.py (Abaqus), mirroring the fairing
+DOE schema (samples[].defect_params + metadata) so downstream tooling is reusable.
+
+Defect types (memo §4):
+  - healthy            — no defect (baseline)
+  - weld_flaw          — porosity / lack-of-fusion on a weld line (located ON a weld)
+  - microcrack         — thermal-cycle micro-cracking (local stiffness loss)
+  - insulation_debond  — foam/MLI separation (same interface type as skin-core debond)
+
+Operating point (memo §2-§3):
+  - temperature_K : 20 K (LH2 cryo) or 293 K (RT), stratified
+  - pressure_MPa  : operational internal pressure
+
+Usage:
+  python src/generate_cryotank_doe.py --n_samples 1 --output doe_cryotank_1sample.json
+  python src/generate_cryotank_doe.py --n_samples 100 --output doe_cryotank_100.json
+  python src/generate_cryotank_doe.py --healthy_only --output doe_cryotank_healthy.json
+"""
+
+import argparse
+import json
+
+import numpy as np
+
+try:
+    from scipy.stats.qmc import LatinHypercube
+    _HAVE_QMC = True
+except Exception:  # scipy missing / too old — fall back to uniform RNG
+    _HAVE_QMC = False
+
+# ---- Geometry bounds (design assumption; see memo §1) ------------------------
+# Al-Li cylindrical LH2 tank, barrel section (representative, dimensions non-public).
+THETA_RANGE = (5.0, 85.0)        # deg, margin from symmetry edges (1/12--1/6 sector)
+Z_RANGE = (400.0, 3600.0)        # mm, along barrel, avoid dome junctions
+WELD_THETA_DEG = [15.0, 45.0, 75.0]  # longitudinal weld-line angular positions
+
+# ---- Defect size tiers (mm): (name, r_min, r_max, fraction) ------------------
+SIZE_TIERS = [
+    ("Small", 20.0, 60.0, 0.35),
+    ("Medium", 60.0, 120.0, 0.40),
+    ("Large", 120.0, 220.0, 0.25),
+]
+
+# ---- Operating point (memo §2-§3; design assumptions) ------------------------
+TEMPERATURE_K = [20.0, 293.0]    # LH2 cryo / room-temperature
+PRESSURE_MPA_RANGE = (0.2, 0.5)  # operational internal pressure
+
+DEFECT_TYPES = ["weld_flaw", "microcrack", "insulation_debond"]
+
+
+def _lhs(n, d, seed):
+    """Latin Hypercube samples in [0,1)^ (n x d); uniform fallback if no scipy.qmc."""
+    if _HAVE_QMC:
+        return LatinHypercube(d=d, seed=seed).random(n)
+    rng = np.random.default_rng(seed)
+    return rng.random((n, d))
+
+
+def _pick_tier(u):
+    """Map u in [0,1) to a size tier by cumulative fraction, return (name, rmin, rmax)."""
+    acc = 0.0
+    for name, rmin, rmax, frac in SIZE_TIERS:
+        acc += frac
+        if u < acc:
+            return name, rmin, rmax
+    name, rmin, rmax, _ = SIZE_TIERS[-1]
+    return name, rmin, rmax
+
+
+def generate(n_samples, seed=42, healthy_only=False, defect_types=None):
+    """Return a DOE dict: {'samples': [...], 'metadata': {...}}."""
+    defect_types = list(defect_types) if defect_types else list(DEFECT_TYPES)
+    samples = []
+
+    if healthy_only:
+        for i in range(n_samples):
+            t_u = _lhs(n_samples, 2, seed)[i]
+            samples.append({
+                "id": i,
+                "job_name": "CryoTank_Healthy_%04d" % i,
+                "defect_params": {
+                    "defect_type": "healthy",
+                    "temperature_K": float(TEMPERATURE_K[int(t_u[0] * len(TEMPERATURE_K)) % len(TEMPERATURE_K)]),
+                    "pressure_MPa": round(float(PRESSURE_MPA_RANGE[0] + t_u[1] *
+                                                (PRESSURE_MPA_RANGE[1] - PRESSURE_MPA_RANGE[0])), 4),
+                },
+            })
+        return {"samples": samples, "metadata": {"n_samples": n_samples,
+                                                 "healthy_only": True, "seed": seed}}
+
+    # 5 design dims: [defect_type, theta, z, size, operating(temp,press packed via 2 dims)]
+    u = _lhs(n_samples, 6, seed)
+    for i in range(n_samples):
+        dtype = defect_types[int(u[i, 0] * len(defect_types)) % len(defect_types)]
+
+        # weld_flaw sits ON a weld line; others use continuous theta
+        if dtype == "weld_flaw":
+            theta = float(WELD_THETA_DEG[int(u[i, 1] * len(WELD_THETA_DEG)) % len(WELD_THETA_DEG)])
+        else:
+            theta = round(float(THETA_RANGE[0] + u[i, 1] * (THETA_RANGE[1] - THETA_RANGE[0])), 2)
+
+        z = round(float(Z_RANGE[0] + u[i, 2] * (Z_RANGE[1] - Z_RANGE[0])), 1)
+        tier, rmin, rmax = _pick_tier(u[i, 3])
+        radius = round(float(rmin + u[i, 4] * (rmax - rmin)), 1)
+        temp = float(TEMPERATURE_K[int(u[i, 5] * len(TEMPERATURE_K)) % len(TEMPERATURE_K)])
+        press = round(float(PRESSURE_MPA_RANGE[0] + (u[i, 5] * 7.0 % 1.0) *
+                            (PRESSURE_MPA_RANGE[1] - PRESSURE_MPA_RANGE[0])), 4)
+
+        samples.append({
+            "id": i,
+            "job_name": "CryoTank_%s_%04d" % (dtype, i),
+            "defect_params": {
+                "defect_type": dtype,
+                "theta_deg": theta,
+                "z_center": z,
+                "radius": radius,
+                "size_tier": tier,
+                "temperature_K": temp,
+                "pressure_MPa": press,
+            },
+        })
+
+    return {"samples": samples,
+            "metadata": {"n_samples": n_samples, "seed": seed,
+                         "defect_types": defect_types,
+                         "size_tiers": [t[0] for t in SIZE_TIERS]}}
+
+
+def main():
+    p = argparse.ArgumentParser(description="Cryogenic hydrogen-tank DOE generator")
+    p.add_argument("--n_samples", type=int, default=1)
+    p.add_argument("--output", type=str, default="doe_cryotank_1sample.json")
+    p.add_argument("--seed", type=int, default=42)
+    p.add_argument("--healthy_only", action="store_true", default=False)
+    p.add_argument("--defect_types", nargs="+", default=None, choices=DEFECT_TYPES)
+    args = p.parse_args()
+
+    doe = generate(args.n_samples, seed=args.seed,
+                   healthy_only=args.healthy_only, defect_types=args.defect_types)
+    with open(args.output, "w") as f:
+        json.dump(doe, f, indent=2)
+    print("Wrote %d samples -> %s" % (len(doe["samples"]), args.output))
+
+
+if __name__ == "__main__":
+    main()

tests/test_cryotank_doe.py

@@ -0,0 +1,42 @@
+"""Smoke tests for the hydrogen-tank DOE generator + FEM scaffold (no Abaqus)."""
+
+import importlib
+
+
+def test_cryotank_doe_generate_one():
+    m = importlib.import_module("generate_cryotank_doe")
+    doe = m.generate(1, seed=42)
+    assert doe["metadata"]["n_samples"] == 1
+    assert len(doe["samples"]) == 1
+    dp = doe["samples"][0]["defect_params"]
+    assert dp["defect_type"] in (["healthy"] + m.DEFECT_TYPES)
+    assert dp["temperature_K"] in m.TEMPERATURE_K
+
+
+def test_cryotank_doe_batch_and_bounds():
+    m = importlib.import_module("generate_cryotank_doe")
+    doe = m.generate(50, seed=7)
+    assert len(doe["samples"]) == 50
+    for s in doe["samples"]:
+        dp = s["defect_params"]
+        assert m.THETA_RANGE[0] <= dp["theta_deg"] <= m.THETA_RANGE[1] or dp["defect_type"] == "weld_flaw"
+        assert m.Z_RANGE[0] <= dp["z_center"] <= m.Z_RANGE[1]
+        assert dp["radius"] > 0
+
+
+def test_cryotank_healthy_only():
+    m = importlib.import_module("generate_cryotank_doe")
+    doe = m.generate(3, healthy_only=True)
+    assert all(s["defect_params"]["defect_type"] == "healthy" for s in doe["samples"])
+
+
+def test_cryotank_dataset_dry_run_plan():
+    m = importlib.import_module("generate_cryotank_dataset")
+    # cryogenic operating point → E should be scaled up vs RT
+    plan = m.resolve_plan({"defect_type": "healthy", "temperature_K": 20.0, "pressure_MPa": 0.3})
+    assert plan["material"]["E_MPa"] > m.MAT_ALLI["rt"]["E_MPa"]
+    assert plan["loading"]["thermal"]["delta_K"] < 0  # cooling from RT reference
+    rt = m.resolve_plan({"defect_type": "weld_flaw", "temperature_K": 293.0,
+                         "pressure_MPa": 0.3, "theta_deg": 45.0, "z_center": 2000, "radius": 80})
+    assert rt["material"]["E_MPa"] == m.MAT_ALLI["rt"]["E_MPa"]
+    assert rt["defect"]["type"] == "weld_flaw"