Orz日記 by Akio Morita

_ [SAD]分離すべき関数

機能毎のSADScript関数の切り離しを行うために汎用で参照される関数を 分離する必要がある

取り合えず、現時点で判明している範囲は次のとおり

• src/tfFileIO_.c
  • cfexptilde_
  • expand_tilde
• src/tfDynamicLink_.c
  • dl_* MI backend
• src/tfBuildInfo_.c.in
  • sadpkgdir_
  • check_build_info
  • sad_pkg_root,build_info_table等のデータベース

作業方針

• cfexptilde_は、src/utils.cに収容(Fortran用だし)
• expand_tildeは、src/sim/sad_api.cへの収容が妥当かな？
• dl_* MI backendは、src/sim/dynl-mi.c辺り？(プロトタイプヘッダは、既にsrc/sim/dynl.hがある)
• BuildInfoに関しては、Feature系同様に独立したframeworkにする(src/buildinfo.[ch]かな？)
• sadpkgdir_は、src/utils.cにてframework呼び出しで実装すべき

多分、src/sim/unix_env_.cの中身も src/utils.cと同じく Fortranコードへのユーティリティー類なのでマージすべきかな？

_ [SAD]SAD->MADX変換

変換精度が出ないエレメントをEmittance[]からtransfer matrixを切り出して matrixで埋め込む作業を始めたのだが、絶賛嵌まり中

SOL内部でSOL/MULT/SOL/MULT/...なチェーンが組んであるとき、 SOLは厚みが無いので、次のMULTの入り口と同じ状態にないと いけないはずだが、Emittance[]から切り出した ClosedOrbitとTransferMatricesは、一つ前のMULTの入り口と同じ値に見える。 多分、SOL内部の計算がtturneからtsoleへ移譲されているが、 演算結果の格納形態が完全互換になっていないように思える。 こいつは、ソースコードを要確認。

取り合えず、次のエレメントを探すコードをハックしてSOL内のMULTを matrixへ展開してみたのだが...結果が全然違う Orz

何が起こっているのだろう？

取り合えず、次の課題は...

• matrixへの展開コードの健全性を別のエレメントで検証する
• SOL内部の転送行列とビーム行列の変換と格納コードを検証する
• 現在省略している DRIFTエレメントのパラメータ展開をサポートする(SOL/DRIFTに転送行列を渡す必要が...)
  • 別案としては、軌道と転送行列を返すハッシュ関数への参照を渡すと言う解も有り

_ [SAD]TableFormバグってる

です。SAD->MADX変換器のデバッグで転送行列を整形して表示するのに 使っていたら、指数部の下一桁が欠落するバグ (例えば、1.23456e-10が、1.23456e-1に化ける) が原因で、思いっきり時間を食った...Orz

つー訳で、amorita branchはRevision 2347で修正入れました。

あと、Emittance[]の返すTransferMatricesとClosedOrbitだが、 どうやらSOL領域は一個ずれているようだ。 原因は、src/tsole.fとsrc/tturne.fでの転送行列(iatr)と軌道(iacod)の 格納法の問題で、src/tturne.fではl番目のエレメントの 入り口の状態をリストのl番目の要素に保存したあとにエレメントの 転送行列を作用させているのに対して、src/tsole.fではl番目の エレメントの転送行列を作用させた結果をリストのl番目の要素に 保存しているので、SOLで囲まれた境界内部では、TransferMatricesや ClosedOrbitのインデックスがズレている。 これと同時に計算されるtwiss配列の更新では、 src/tsole.fではリストのl+1番目の要素へ格納するので整合する。 多分、格納コードをコピー＆ペーストしたが 保存と転送の順番を入れ替えに対応してコードを書き換えるのを 忘れたと推測される。

修正は、Revision 2348で修正をコミットした。

アーキテクチャ的には、SOL領域内部の写像は ソレノイド磁場が0となる極限ではSOL領域外と一致するので、 単一のコードで処理すべきだと思う。

_ [SAD]SAD->MADX変換

昨日の問題だが、原因が判明

• Emittance[]から得た転送行列をmatrixエレメントへ変換する際に、ROTATEパラメータによる不要な座標回転変換を加えていた
• SOL内部のインデックスのズレにより、Emittance[]から得た転送行列と軌道との関係が、正しくないのでアフィン変換の定数項が正しく計算出来ていなかった
• SOL内部のインデックスのズレに対する、Quick Hackの際にSOL境界を越えてNext Elementを探索したので、FBLX(何もしないMULT)を正しく単位行列に変換出来なかった

一つはこちら側のミスだが、残りはSADの実装バグに起因する Orz

Emittance[]で取り出した転送行列でエレメントを置換するルーチンは完成

SOL系のエレメントを転送行列近似に切り替えて様子を見ると

• HERでは、分散関数と軌道の近似が劣化
• LERでは、分散関数と軌道の近似が改善(主に垂直方向)

という、悩ましい結果に...

QCSL側で、まっすぐ通過している(LER) or 傾いている(HER)の違いが 現れているのだろうか？

別の可能性としては、Emittance[]の返す転送行列が先頭から 各エレメントまでのものなので、エレメント毎の行列を 行列積で計算することによる誤差混入とかSAD-MADX間の 座標系定義の差異とかが疑わしい。

座標系定義に関してはドキュメントをあたって定義を突き合わせてみる。

全部のエレメントを転送行列で展開しても軌道に残差が残る Orz

転送行列とCODをアフィン変換へ展開して、MADX内部で再度 転送行列とCODを計算する過程の演算精度が原因か？

次のステップは、SAD上でアフィン変換展開して後でCODの再構築での 再現精度を検証することかな？

_ [SAD]SAD->MADX変換

転送行列とCODからのアフィン変換導出と、アフィン変換列からの転送行列と CODの再計算は、それなりの精度で出来ることを確認した。 残る問題は、浮動小数点の text dump時の精度落ちと MADXの内部計算での誤差の蓄積と言うことになる。

• text dump時の精度落ち
  • Write[]/Get[]して調べたが、問題となるほど悪化していない

_ [MADX]3.04.40のkickerの tiltの向きが違う件はバグ認定されたらしい

最初は、「バグではない仕様です」と言っていたらしいが 流石にモジュール間で軌道が一致しないのことに気がついたのか バグ認定したらしい。(バグ追跡システムが無いので人を介しての伝聞形)

現象としては、tiltを付けたkickerを含むビームラインでの 軌道のtilt依存性が

• MADX-3.04.40のtwissとマニュアルの記述が逆符号になる
• MADX-3.04.40のtwissとtrackで逆符号になる
  • つまり、trackの結果は正常
• MADX-3.04.40とMADX-3.03.25間でtwissの結果の符号が異なる

状況証拠的には、内部実装の変更時にエンバグしたというシナリオが 一番しっくりくるのだが、最初に「kickerのtiltの符号が 間違ってない？」と報告した時は、「kickerは特殊なエレメントだから これで正しい」との返事が返ってきたらしい。

開発現場では、何をテストして正常だと言ったのだろう...Orz

つー訳で、このバグが直るとSAD->MADX変換のworkaroundを外さないといけない...

_ [SAD]SAD->MADX変換

全てのSADエレメントをmatrixに変換したモデルで軌道が再現しない件だが、 SOL領域のエレメントに僅かな隙間が開いてそこへドリフトが 挿入されることを発見した。

原因は、当該領域ではLINE["S", pos + 1]が LINE["S", pos] + LINE["L", pos]と等しく無いため。 多分、SOL領域ではトラッキングはSOLの軸上で行われるが 軌道は一般には斜めに走るので、軸上で見た距離(エレメント長Lの積算)と 軌道長の積算(S)に差が生じているためと思われる。

エレメント分割管理と、エレメントの位置情報の計算を2重台帳 (SAD上のSとエレメント軸上の距離の積算)化して対応したが、 なんかたちの悪いバグが残っているようなので、調査中...

LERの変換結果が途中からおかしくなってる Orz

富士のBSWFRP以降の軌道に分散関数のパターンが見える... アフィン変換自体はBSWFLPと比べても異常が無い。 あれこれ調べた結果、Emittance[]が出す転送行列は 6x6で、 軌道も加速込みの平衡軌道になっているが、 初期値に使ったのがTwiss[]の結果なので、 SOLをくぐる所でDZがずれた結果、富士の加速空洞で加速が起こり 加速空洞の次のBENDである BSWFRPから分散関数の分だけ 軌道のズレを生じている...

エレメントの内点での計算を考えると、6x6なTwiss[]が欲しいよね

_ [SAD]コードの整理作業

下記のファイル中からファイルスコープ外から参照される関数を撤去した。

• src/tfFileIO_.c
• src/tfDynamicLink_.c
• src/tfBuildInfo_.c.in

BuildInfoに関しては、src/buildinfo.hにある公開APIから検索可能にした。

_ [SAD]Cygwin環境

久々に、Cygwin環境でコンパイルしてみた。 で、色々問題を見つけたのでまとめてみる。

• implicit declaration of function 'strdup'(src/random_driver.c,src/feature.c,src/sim/sad_api.c)
  • string.hをincludeしてるはずなのに、なぜかでる

• 'M_PI' undeclared (first use in this function)(src/random_plugin_sad.c)
  • math.hをincludeしてるはずなのに、なぜかでる

• implicit declaration of function 'hypot'(src/sim/fortran_math_hypot_.c)
  • ISO C99の数学関数(math.h)

• implicit declaration of function 'asinh'(src/sim/fortran_math_asinh_.c)
  • ISO C99の数学関数(math.h)

• implicit declaration of function 'acosh'(src/sim/fortran_math_acosh_.c)
  • ISO C99の数学関数(math.h)

• implicit declaration of function 'atanh'(src/sim/fortran_math_atanh_.c)
  • ISO C99の数学関数(math.h)

• implicit declaration of function 'posix_openpt/grantpt/ptsname'(src/sim/fortran_io_.c)
  • IEEE Std 1003.1-2001(POSIX.1)関数

• ucontext.hが無い(src/sim/sad_signal.h)
  • SUSv2由来

• implicit declaration of function 'vsnprintf'(src/sim/sad_api.c,src/tfProcess_.c)
  • ISO C99関数

• implicit declaration of function 'snprintf'(src/tfNetworkIO_.c,src/tfProcess_.c)
  • ISO C99関数

• implicit declaration of function 'isnan'(src/utils.c)
  • ISO C99関数

• implicit declaration of function 'tempnam'(src/utils.c)
  • ISO C90関数

• implicit declaration of function 'mkstemp'(src/sim/fortran_io_.c,src/tfFileIO_.c)
  • IEEE Std 1003.1-2001(POSIX.1)関数

• implicit declaration of function 'realpath'(src/tfFileIO_.c)
  • IEEE Std 1003.1-2001(POSIX.1)関数

• implicit declaration of function 'fmin/fmax'(src/tfDateTAI_.c)
  • ISO C99の数学関数(math.h)

• implicit declaration of function 'crypt'(src/tfCrypt_.c)
  • IEEE Std 1003.1-2001(POSIX.1)関数

• implicit declaration of function 'putenv/unsetenv'(src/Process_.c)
  • IEEE Std 1003.1-2001(POSIX.1)関数

とにかく、ISO Cと POSIXまわりで足りないものがあるような...

_ [LHC]運休中 Orz

冷凍系の作業で水曜日まで運休中

_ [SAD]newlibだめだめ Orz

Cygwinのエラーの件だが、 GCCのソースを読み進めて行くと、-std=オプションを使って ISO Cを指定したり(例えば-std=c99)、-ansiオプションが 設定された場合、cppにてマクロ__STRICT_ANSI__が 自動的に定義される(参照: gcc/c-cppbuiltin.c)。

ここまでは問題無いのですが、Cygwinの newlib付属のヘッダー群は、 __STRICT_ANSI__定義時にそれなりの数の関数の定義が無効になる と言う動作をしています。

FreeBSDのヘッダー(newlibの元ネタは BSD libc)を覗いてみると /usr/include/sys/cdefs.h当たりで !defined(__STRICT_ANSI__) || __STDC_VERSION__ >= 199901 といった使い方をしています。 意味的には、__STDC_VERSION__が 199901未満 (つまり、ISO C99以前)かつ、 __STRICT_ANSI__が定義されているときに、 以降のヘッダー宣言が無効になるというものです。

__STRICT_ANSI__定義時は、厳格なANSI規格に適合しない コードを除外しているわけだが、ISO C99な関数の宣言まで無効になっている ということは、newlibヘッダーは最新の規格に追従出来てないということ。

と言う訳で、CygwinでコンパイルできないのはSAD側のソースの問題ではなく Cygwinが標準規格に追従出来ていないのが主原因と言える。

_ [SAD]閏秒データベース更新

現在、SADでの UTCと SAD epoch timeの相互変換には使っている libtai用の 閏秒データベースに 2008-12-31 23:59:60 UTCを挿入した。

これは、 閏秒対応の unix timeで運用している計算機のみが影響を受ける話だが、 年末までにはバックポートを済ますべきかな...

_ [SAD]SAD->MADX変換

変換時のモデリングと軌道・分散関数の再現性だが、 もっとも高い再現性を得る手法はEmittance[]から得られる TransferMatrices/ClosedOrbitに基づいてアフィン変換列を 作り、それを直接MADXのmatrix要素に使う手法である。

それを基準に、各要素をMADXの要素に置換した場合の 軌道・分散関数の再現性への影響

• 軌道再現性
  • DRIFT/QUAD/SEXTは大きな影響無し
    • ただし、LERの垂直シケインで軌道エラーが見える(matrix展開でも若干現れるので、MADXとSADの軌道計算モデルの差異が原因か？)
  • BEND/MULT等のK0パラメータを含むもの(hkickerへ展開する)は、NG
  • SOLの境界座標変換は影響無し
  • SOL/DRIFTをsolenoidへ変換すると再現性が桁で悪化する
  • SOL/MULTをsolenoid+multipole+hkickerへ変換するとSOL/DRIFTより影響が大きい(2〜3倍？)
    • 垂直方向の再現性への影響はSOL/MULT >> SOL/DRIFT〜SOL

• 分散関数再現性
  • BEND/MULT等のK0パラメータを含むもの(hkickerへ展開する)は、NG
  • DRIFT/QUAD/SEXTも地味に影響が出る＆影響の出方がバラバラ(分散フリーな直線部でμmオーダー)
  • SOLの境界座標変換は影響無し
  • SOL/DRIFTをsolenoidへ変換すると再現性が桁で悪化する
  • SOL/MULTをsolenoid+multipole+hkickerへ変換するとSOL/DRIFTより影響が大きい(2〜3倍？)

SOLの境界座標変換に関しては、Emittance[]の計算エンジンsrc/tsconv.f から持ってきているので当然の結果と言える。

現時点では、SOL内部のエレメント(DRIFT/BEND/QUAD/MULT)に関しては SAD側のモデルの転送行列を移植する以外の再現手段が無い。 (物理的な内容を反映したモデリングでは、軌道すら再現できない)

BEND/MULTで K0を含む場合の分散関数の発生や軌道の偏向の再現性に難がある。 おそらく、MADX内部の kickerのモデリングの問題がある＆ SADでは ANGLE依存の出入り口の座標変換と ANGLE+K0依存の軌道偏向計算を 取り扱っているためと思われる。

_ [SAD]64bitSADの可能性

gfortran 4.4でI/OのencodingにUTF-8がサポートされ、内部の文字列を UCS-4(kind = 4)で扱うことが出きるようになったらしい。

原理的には、8倍精度実数(kind = 32)と128bit整数(kind=16)が有れば、 real = 2 * integer = 8 * characterの関係を維持したままLP128化出来る。

問題は、8倍精度実数なんてものを実装したハードウェアが実在しない辺りかな...

4倍精度実数関しては、IEEE754rのドラフトベースのハードウェア実装が SunのSPARCやIBMのPOWERに有るらしいので、gfortranを改造して UCS-2(kind = 2)な内部表現をサポートしてもらえば、 4倍精度実数(kind = 16)と64bit整数(kind=8)と組み合わせて real = 2 * integer = 8 * characterの関係を維持したままLP64化出来る。

_ [Admin][tDiary]ツッコミSPAM

頑張って投稿する連中(たぶんBOT)が多いので、色々対策する必要を 迫られる訳だが、なぜが 2008-02-20の記事が集中的に狙われる。

不思議に思って調べてみたのだが、検索エンジンからの紐が 一番多く繋がっている記事を狙うアルゴリズムのようだ...

_ [SAD]SAD<->MADX変換

SAD->MADX変換は、これ以上進捗しそうに無いので、MADX->SAD変換も 含めて手元にコードをSAD Source Archiveへ置いてきた。

KEKB latticeの変換に関しては、SADとMADXのモデルが違いすぎるので MADX側を改造しない限り再現は無理な気がする。 現時点で把握している、再現上の問題は...

• 大きなoffsetを持った軌道(IR部のsolenoid領域)
  • 近軸近似の精度？(MADXのTwiss.Fは2次項まで使っている)
• K0付きBEND/MULT
  • MADXのrbend/multipoleにK0(設計軌道を曲げずに軌道を曲げる)パラメータが無い
  • kicker要素は、dipole kickのみで edge focusを含まない

さらに、MAP要素を含めるとSAD->MADX変換システムを構築することは、 次の命題と同値と思われる。

チューリング完全な言語で記述される任意の写像関数を必要な精度で近似するMADX sequenceを組み立てる

MAP要素ではSADScriptが、一般の要素ではFortranコードが入力 とするある種のコンパイラを作る話になる。

最低でも、MADX側の element typeは任意のsequenceを組み立てたときに $R^6 \rightarrow R^6$写像関数の空間を埋め尽くせるだけの多様性が必要。

_ [SAD][MADX]SAD->MADX変換

「MADX-PTCなら全てのエレメントがきちんと実装されてる」なんて言われたんで、 MADX-PTC試してるんですけど...

++++++ warning: ptc_twiss:  DA got unstable after normal
++++++ Error: seterrorflag : Errorcode: 10   Reported from ptc_twiss :
++++++ Error: seterrorflag : Description: DA got unstable after normal

とかいって、全然計算してくれません

twissが計算できることを計算できないptc_twissって使えね〜

_ [SAD][MADX]SAD->MADXP変換の問題

色々有るけど、現時点で致命的なのは

• matrix要素が使えない
  • QUADのlinear fringeが変換できない
  • thinなBENDが変換できない
  • SOLの境界での座標変換が記述できない
  • DX,DYパラメータが記述できない

辺りかな...

PTCはどんな加速器でも記述できると言う人がいるが、 Front-end processorにMADXを使う限り MADXよりも 記述能力が低い罠 Orz

追記(2008-09-24 16:00)

Frank Schmidtに、 「matrixは擬似エレメントだ、PTCに対応するエレメントは無い」と 言われた。と言う訳で、SAD->MADXP変換は終了〜〜〜

_ [SAD][YaSAI]省略可能な構文要素の記述法

文法上、構文要素が省略可能かつ省略時に暗黙のうちに仮定される 構文要素がある場合、単純にOptional記述で文法を書き下すと PEGパーサーが組み立てる抽象構文木(AST)は、省略が有った場合と 露に記述された場合とで異なるものになる。 特に再帰の無い文法要素(例えば、SADの CompoundExpression)では、 AST評価器で逐一リーフノードの型を検査する必要がある。

具体的な例としては、

Compound <- (Expr? OpCompound)+ Expr?

な文法では、CompoundのASTの枝にはExprと OpCompoundが順不同に現れる。 ここで、Exprが省略された時に、そこに暗黙のExprが 仮定される場合、次の文法を使うことでAST評価器を単純化できる。

Compound <- ((Expr / Null) OpCompound)+ (Expr / Null)
Null     <- &.?

ポイントはNullの定義で、これは入力を消費せずに 常に真となるルールである。つまり、Exprがマッチしなければ ASTへNullノードを挿入して次のルールの評価に進む。 完成したCompoundのASTでは、OpCompoundとOpCompoundの間には、 常にExprかNullが挟まれるので、 OpCompoundであるかどうかの検査が不要になる。

トレードオフは、ASTがNullノードを含む分肥大化することと、 構文解析器が余計な状態を持つことかな？

メリットとしては、AST評価器でリーフノードをたどるループが 単純化することと、省略時に仮定する暗黙の構文要素の種別を 文法に露に表現できること。

取り合えず、CompoundExpressionと Consルールで動作を実証したので、 こいつを使って Alternativesとかでのオペランド省略時に仮定される 暗黙のNullシンボルをサポートする予定。

_ [雑記]TODOタスクリスト

• LHC関連
  • Crab Cavity Frequency Sweep Simulationの実装とベンチマーク
  • Higher Order Multipole付きのCrab Cavity Ramping Simulation(パラメータ待ち)
  • Crab Crossing時の Injection Mismatch Simulation
    • 想定されるInjection Mismatchの関数形の設定
    • 評価関数の設定
• SAD関連
  • SAD->MADX変換エンジン開発(凍結)
  • Tcl/Tk vs xproto-7.0.13問題対応
  • Parallel Algorism Wrapperのテンプレートクラスの改良
  • YaSAIまわり
    • cstdatomicを使ったReader/Writer Lockの実装と Multi-thread Safe化(最初はSymbolIDクラスから)
    • mapをunordered_mapへ置き換えたベンチマークとコード書き換え