dynamicsoar's log

主に研究関係のメモ

Fluent UDF: NODE_POS_UPDATED で付与された MOVED フラグは NODE_POS_MOVABLE コマンドでリセットできる

問題編

DEFINE_GRID_MOTION を使って板とかを曲げたりしたいとき、メッシュ内の各 node は複数の face にまたがっている可能性があるので、UDF マニュアルのサンプルにあるように、

#include "udf.h"

DEFINE_GRID_MOTION(my_dyn_mesh, domain, dt, time, dtime)
{
  Thread *tp = DT_THREAD(dt); // tp = thread pointer
  face_t f; // face index
  Node *np; // np = node pointer

  begin_f_loop(f,tp)
  {
    f_node_loop(f,tp,n)
    {
      np = F_NODE(f,tp,n);

      if ( NODE_POS_NEED_UPDATE(np) )
      {
        NODE_POS_UPDATED(np); // Update the flag
        
        do_something;
      }
    }
  }
  end_f_loop(f,tp);
}

というように、

  1. if( NODE_POS_NEED_UPDATE(np) ) で「既に訪れたかどうか」のフラグを調べて、
  2. NODE_POS_UPDATED(np); コマンドで「訪れた」という情報に書き換える。

しかし、「計算開始直後など、時間進行を伴わずに全ノードをループして何かをしたい」ということがある(自分はあった)。このとき、何も考えずにループ構造をコピペしてしまうと、1回目のループで既にフラグが「既に訪れた」ことになっているので、2回目のループでは「全ての」node をスキップしてしまう。つまり何も起きない。この問題にしばらく悩んでいた。つまり、

DEFINE_GRID_MOTION(my_dyn_mesh, domain, dt, time, dtime)
{
  Thread *tp = DT_THREAD(dt); // tp = thread pointer
  face_t f; // face index
  Node *np; // np = node pointer

  // First loop
  begin_f_loop(f,tp)
  {
    f_node_loop(f,tp,n)
    {
      np = F_NODE(f,tp,n);

      if ( NODE_POS_NEED_UPDATE(np) )
      {
        NODE_POS_UPDATED(np);
        
        do_something;
      }
    }
  }
  end_f_loop(f,tp);

  // Another loop
  begin_f_loop(f,tp)
  {
    f_node_loop(f,tp,n)
    {
      np = F_NODE(f,tp,n);

      if ( NODE_POS_NEED_UPDATE(np) )
      {
        NODE_POS_UPDATED(np);
        
        do_something_else; // This doesn't work!
      }
    }
  }
  end_f_loop(f,tp);
}

ということだ。…と思っていた。調査編でわかるように、実際には単に「訪れたかどうか」ではなかったのだが…。

調査編

まず、フラグの実体はなんなのか。 適当なフォルダ(C:\Program Files\ANSYS Inc\v191\fluent\fluent19.1.0\ とか)において、Devas で NODE_POS_NEED_UPDATENODE_POS_UPDATED を検索すると、dynamesh_tools.h というヘッダファイルにヒットする。

/* node marks */
#define NODE_NULL       0x0000
#define NODE_FIXED      0x0001
#define NODE_MOVABLE    0x0002
#define NODE_CONTACT    0x0004
#define NODE_MOVED      0x0008
#define NODE_ADAPT      0x0010
#define NODE_DEFORMED   0x0020
<中略>
#define NODE_POS_MOVABLE(v) \
  (NODE_CLEAR_MARK(v, NODE_MOVED | NODE_FIXED | NODE_ADAPT), \
   NODE_SET_MARK(v, NODE_MOVABLE))
#define NODE_POS_NEED_UPDATE(v) (NODE_MARK(v) & NODE_MOVABLE)
#define NODE_POS_UPDATED(v) \
  (NODE_CLEAR_MARK(v, NODE_MOVABLE | NODE_FIXED | NODE_ADAPT), \
   NODE_SET_MARK(v, NODE_MOVED))
#define NODE_POS_UPDATED_P(v) (NODE_MARK(v) & NODE_MOVED)
<略>

とある。 まず、フラグは TURE/FALSE のような2値ではなく、複数の値をとる。それがソース冒頭の NODE_NULL, NODE_FIXED, ... というやつだ。 そしてこのフラグに関連するコマンドは「現在のフラグ状況を調べるコマンド」と「フラグをセットするコマンド」の2種類に分かれる。 ここでは MOVABLE と MOVED 関連のみを示したが、実際には他のフラグに対応するコマンドもある。

現在のフラグ状況を調べるコマンド: - NODE_POS_NEED_UPDATE: MOVABLE かどうかを調べる - NODE_POS_UPDATED_P: MOVED かどうかを調べる

フラグをセットするコマンド: - NODE_POS_MOVABLE: MOVED 等のフラグを CLEAR_MARK して、MOVABLE をセットする - NODE_POS_UPDATED: MOVABLE 等のフラグを CLEAR_MARK して、MOVED をセットする

となっている。調べるコマンドは基本的に suffix として _P が付くのだが、なぜか MOVABLE だけは MOVABLE_P ではなく NEED_UPDATE となっていて一貫性に欠ける。

解決編

したがって、解決策としては、最初のループの後でフラグリセット用のループを回して、そこで NODE_POS_MOVABLE を実行する、ということになる。つまり、

DEFINE_GRID_MOTION(my_dyn_mesh, domain, dt, time, dtime)
{
  Thread *tp = DT_THREAD(dt); // tp = thread pointer
  face_t f; // face index
  Node *np; // np = node pointer

  // First loop
  begin_f_loop(f,tp)
  {
    f_node_loop(f,tp,n)
    {
      np = F_NODE(f,tp,n);

      if ( NODE_POS_NEED_UPDATE(np) ) // Check the MOVED flag
      {
        NODE_POS_UPDATED(np); // Set MOVED flag, clear any other flags
        
        do_something;
      }
    }
  }
  end_f_loop(f,tp);

  // Loop for resetting the MOVED flag
  begin_f_loop(f,tp)
  {
    f_node_loop(f,tp,n)
    {
      np = F_NODE(f,tp,n);
      NODE_POS_MOVABLE(np); // Clear MOVED flag (and other flags if any), set the MOVABLE flag
    }
  }
  end_f_loop(f,tp);

  // Another loop
  begin_f_loop(f,tp)
  {
    f_node_loop(f,tp,n)
    {
      np = F_NODE(f,tp,n);

      if ( NODE_POS_NEED_UPDATE(np) ) // Check the MOVED flag
      {
        NODE_POS_UPDATED(np); // Set MOVED flag, clear any other flags
        
        do_something_else; // This works!
      }
    }
  }
  end_f_loop(f,tp);
}

ということだ。

ところで、こうやってループを複数書いてると、end_f_loop(f,tp); を忘れてハマりがちなので注意したい(ハマった)。

ANSYS Meshing と Fluent (or Fluent Meshing) の表示する maximum skewness や minimum orthogonal quality が微妙に違う(原因不明)

再現手順

  1. ANSYS Meshing でメッシュを作成し、maximum skewness を確認する。今回は 0.944 くらいだった。まぁ < 0.95 だし、いけるかな、と。
  2. Fluent 用に .msh を export する(Workbench で Update とするのと同じことだが、自分はこちらを好む)。
  3. Fluent や Fluent Meshing(Fluent launcher で Meshing mode を ON)で読み込む。
  4. TUI の /mesh/check-verbosity を 2 以上にする(デフォルトは 0)。
  5. Mesh quality をチェックする (TUI だと mesh/quality)。すると Maximum Cell Equivolumen Skewness = 0.967 とか出てくる。なんだこれ。minimum orthogonal quality も違う。誤差ってレベルじゃねーぞ

考えられる原因

おそらく

  1. 定義が微妙に違う
  2. 実は .msh に export するときに微妙に形状が変わっている

のどちらかだろう。後者だとさすがにひどいから前者だと思いたいが、もしそうだとしてもなんでだよ感がすごい。もともとは別の会社のアプリだったとかこっちは知ったことじゃない。統合して売ってるんだから整合性とってくれよ…

Matlab で複数のスペースで区切られたデータ (separated by more than one white space) を読み込む → textscan → と思ったが importdata でいいか

Reading in ascii files with white space as delimiter. - MATLAB Answers - MATLAB Central にあった。というかもちろん公式ヘルプもある→ Read formatted data from text file or string - MATLAB textscan - MathWorks United Kingdom。けど、これの存在を知ってないとたどり着けないわけで…

たとえば

f n NODE_X(v) NODE_Y(v) NODE_Z(v)
0 0  -1.5000000e-02   6.1000001e-02   9.4994903e-11
0 1  -1.5000000e-02   6.0000000e-02   2.0000000e-03
0 2  -1.5000000e-02   5.9000001e-02   9.4994903e-11
1 0  -1.5000000e-02   1.0700000e-01   9.4994903e-11
1 1  -1.5000000e-02   1.0800000e-01  -2.0000000e-03
1 2  -1.5000000e-02   1.0900000e-01   9.4994903e-11

みたいなデータを読みたいとする。数値の数値の間に2つ以上のスペースがあるので、ふつうに dlmread とかするとそのスペース自体を 0 とみなされる(っぽい)。これは困る。そこで、

fid = fopen('input.dat')
matrix = textscan(fid, '', 'HeaderLines',1, 'CollectOutput',true);

とすれば行ける。HeaderLines のところの 1 はヘッダーが1行だから。

…と思っていたんだけど、importdata でいいわ。

git completion (for bash) & prompt を入れる

git/contrib/completion at master · git/git · GitHub

ここにある。やるべきことはスクリプト内に全部書いてあるが、一応自分用にメモ。

  1. git-completion.bashgit-prompt.sh をそれぞれ .git-completion.bash.git-prompt.sh の名前で ~/ に保存する。CRLF だと(少なくとも Cygwin terminal だと)動かないので LF にすること。
  2. ~/.bashrc に以下を追加:
# For git-completion & prompt
source ~/.git-completion.bash
source ~/.git-prompt.sh
PS1='[\u@\h \W$(__git_ps1 " (%s)")]\$ '

ANSYS Fluent の autosave 間隔 n は「計算開始(再開)時点から n ステップごと」ではなく「現在ステップ/n が割り切れたら出力」の仕様っぽい

たとえば、非定常計算で、

;; Autosave Every (time step)
file/auto-save/data-frequency
10
solve/set/time-step
0.02
solve/dual-time-iterate
50
30

として 1 s まで計算したとする。こうすると、タイムステップが 10, 20, 30, 40, 50 のときに autosave で .dat ファイルが出力される。計算が終わるタイムステップはもちろん 50 だ。

次に、この計算に続けて(あるいはリスタートして)、

;; Autosave Every (time step)
file/auto-save/data-frequency
100
solve/set/time-step
0.002
solve/dual-time-iterate
2000
30

のように、dt を 1/10 にして、autosave 間隔も 10 倍の 100 に変えるとどうなるか。「計算開始(再開)時点のタイムステップから始まって 100 ステップごと」に autosave されるので、

  • 150
  • 250
  • 350
  • 450
  • 550
  • 650 ...

となると思っていた。違った。実際には

  • 100
  • 200
  • 300
  • 400
  • 500
  • 600 ...

が出力される*1。つまり、「計算開始(再開)時点のタイムステップから始まって 100 ステップごとに保存」されるというのは間違いで、「単に、各時点でのタイムステップを autosave frequency で割って、割り切れたら保存」としているようだ。

やっていないが、おそらく定常計算での iteration 間隔での保存でも全く同様なのだろう。

この仕様はいまの自分の計算的には不便で、残念だ。しかしそれとは別の話として、このソフトにはこういう「やってみるまで挙動がわからない(マニュアルにも書いていないか、書いていても絶望的に見つけづらい)仕様」が多すぎて、つらいところがある。

*1:と思われる。自分がやったのはこの数字ではなかったが、とにかく「現在のステップから」ではなく「割り切れる値」のときになった。

Office Lens で OneDrive へのアップロードが不能に → Office Lens に加えて OneNote, SharePoint からも sign out したら直った

タイトルの通り。

具体的な症状としては、My Files 内の画像アイコンが "Unable to Upload" となっていた。Retry ボタンを押すとすぐに大学の Office 365 の sign in 画面へ飛ばされ、正しい情報を入力しても Office Lens に戻されて、何も起きない、という状況だった。

少しググったところ、「OneDrive のカメラ機能を使いましょう」というMS公式のページがヒットして、「おいふざけんな… Office Lens の意義はなんだよ…」ってなったが、さらにもう少しググると、タイトルの方法(Office Lens だけでなく他の MS アプリからも sign out する)が有効という報告を見かけたのでやってみた。効いた。

なお iOS 標準の Mail アプリで Exchange のメール受信をしていて、これは sign out しなかったが問題なかった。

2019-08-08 追記

また同じ症状になった。このエントリのことを忘れていたのでまたググってしまったのだが、

Office Lens - won't save to O365 OneDrive - Microsoft Community

にあるように、一度アプリを全部閉じてから、iOS の Settings > OneDrive から Clear Account settings を ON にして、再度 Office Lens を開き直して sign in したら解決した。