粒子法による流体解析
termat:http://termat.sakura.ne.jp/
*
* MPS法（moving particle semi-implicit method、粒子法の一種） による流体解析の習作
*
* 参考文献：粒子法(計算力学レクチャーシリーズ) 越塚 誠一 著 丸善出版
* url=http://pub.maruzen.co.jp/book_magazine/book_data/search/4621075225.html
*
* 1.注意事項
* 計算にとても時間がかかります。
* OS=Windows7、CPU=Core i7、 Memory=4GBの環境で、計算に、
* 問題No1=80秒、問題No2=247秒、 問題No3＝335秒を要しました。
*
* 2.操作
* 1)チェックボックスで問題を選ぶと、計算を開始します。
* 2)実時間で2秒分の計算を終えると、計算結果を表示します。
* 3)計算が終了した段階で画面をクリックすると、計算結果を再表示します。
*
* 3.備考
* 計算条件ファイル（CSVファイル形式）は、以下の書式でデータを保存している。
* 1行目-------平均粒子間距離
* 2行目以降----格子タイプ,x座標,y座標
*
/*
* termat:http://termat.sakura.ne.jp/
*
* MPS法（moving particle semi-implicit method、粒子法の一種） による流体解析の習作
*
* 参考文献：粒子法(計算力学レクチャーシリーズ) 越塚 誠一 著 丸善出版
* url=http://pub.maruzen.co.jp/book_magazine/book_data/search/4621075225.html
*
* 1.注意事項
* 計算にとても時間がかかります。
* OS=Windows7、CPU=Core i7、 Memory=4GBの環境で、計算に、
* 問題No1=80秒、問題No2=247秒、 問題No3＝335秒を要しました。
*
* 2.操作
* 1)チェックボックスで問題を選ぶと、計算を開始します。
* 2)実時間で2秒分の計算を終えると、計算結果を表示します。
* 3)計算が終了した段階で画面をクリックすると、計算結果を再表示します。
*
* 3.備考
* データ：格子タイプ,x座標,y座標
*
* 4.更新履歴
* 2014/12/14 メンテナンスを実施
*/
package
{
import flash.display.Bitmap;
import flash.display.BitmapData;
import flash.display.MovieClip;
import flash.display.Sprite;
import flash.events.Event;
import flash.events.MouseEvent;
import flash.net.URLLoader;
import flash.net.URLRequest;
import flash.system.Security;
import flash.text.TextField;
import flash.text.TextFormat;
import com.bit101.components.*;
[SWF(width = "480", height = "480", backgroundColor = "0x000000", fps = "60")]
public class MPS3 extends Sprite {
private var loader:URLLoader;
private var cont:MPSController;
private var scale:Number;
private var colors:Array = [0x0000ff, 0x555555, 0xaaaaaa];
private var radius:Number;
private var time:TextField;
private var bitmap:BitmapData;
private var mc:MovieClip;
private var mc2:MovieClip;
private var ball:Vector.<Metaball>;
private var rec:Vector.<Array>;
private var times:Vector.<Number>;
private var progress:ProgressBar;
private var endTime:Number = 2.0;
private var isThinking:Boolean = true;
private var counter:int = 0;
private var nextTime:Number = 0.02;
private var isRendering:Boolean = true;
private var prob01:String="0,0.012,0.012/0,0.012,0.028/0,0.012,0.044/0,0.012,0.06/0,0.012,0.076/0,0.012,0.092/0,0.012,0.108/0,0.012,0.124/0,0.012,0.14/0,0.012,0.156/0,0.012,0.172/0,0.012,0.188/0,0.012,0.204/0,0.012,0.22/0,0.012,0.236/0,0.012,0.252/0,0.012,0.268/0,0.012,0.284/0,0.028,0.012/0,0.028,0.028/0,0.028,0.044/0,0.028,0.06/0,0.028,0.076/0,0.028,0.092/0,0.028,0.108/0,0.028,0.124/0,0.028,0.14/0,0.028,0.156/0,0.028,0.172/0,0.028,0.188/0,0.028,0.204/0,0.028,0.22/0,0.028,0.236/0,0.028,0.252/0,0.028,0.268/0,0.028,0.284/0,0.044,0.012/0,0.044,0.028/0,0.044,0.044/0,0.044,0.06/0,0.044,0.076/0,0.044,0.092/0,0.044,0.108/0,0.044,0.124/0,0.044,0.14/0,0.044,0.156/0,0.044,0.172/0,0.044,0.188/0,0.044,0.204/0,0.044,0.22/0,0.044,0.236/0,0.044,0.252/0,0.044,0.268/0,0.044,0.284/0,0.06,0.012/0,0.06,0.028/0,0.06,0.044/0,0.06,0.06/0,0.06,0.076/0,0.06,0.092/0,0.06,0.108/0,0.06,0.124/0,0.06,0.14/0,0.06,0.156/0,0.06,0.172/0,0.06,0.188/0,0.06,0.204/0,0.06,0.22/0,0.06,0.236/0,0.06,0.252/0,0.06,0.268/0,0.06,0.284/0,0.076,0.012/0,0.076,0.028/0,0.076,0.044/0,0.076,0.06/0,0.076,0.076/0,0.076,0.092/0,0.076,0.108/0,0.076,0.124/0,0.076,0.14/0,0.076,0.156/0,0.076,0.172/0,0.076,0.188/0,0.076,0.204/0,0.076,0.22/0,0.076,0.236/0,0.076,0.252/0,0.076,0.268/0,0.076,0.284/0,0.092,0.012/0,0.092,0.028/0,0.092,0.044/0,0.092,0.06/0,0.092,0.076/0,0.092,0.092/0,0.092,0.108/0,0.092,0.124/0,0.092,0.14/0,0.092,0.156/0,0.092,0.172/0,0.092,0.188/0,0.092,0.204/0,0.092,0.22/0,0.092,0.236/0,0.092,0.252/0,0.092,0.268/0,0.092,0.284/0,0.108,0.012/0,0.108,0.028/0,0.108,0.044/0,0.108,0.06/0,0.108,0.076/0,0.108,0.092/0,0.108,0.108/0,0.108,0.124/0,0.108,0.14/0,0.108,0.156/0,0.108,0.172/0,0.108,0.188/0,0.108,0.204/0,0.108,0.22/0,0.108,0.236/0,0.108,0.252/0,0.108,0.268/0,0.108,0.284/0,0.124,0.012/0,0.124,0.028/0,0.124,0.044/0,0.124,0.06/0,0.124,0.076/0,0.124,0.092/0,0.124,0.108/0,0.124,0.124/0,0.124,0.14/0,0.124,0.156/0,0.124,0.172/0,0.124,0.188/0,0.124,0.204/0,0.124,0.22/0,0.124,0.236/0,0.124,0.252/0,0.124,0.268/0,0.124,0.284/0,0.14,0.012/0,0.14,0.028/0,0.14,0.044/0,0.14,0.06/0,0.14,0.076/0,0.14,0.092/0,0.14,0.108/0,0.14,0.124/0,0.14,0.14/0,0.14,0.156/0,0.14,0.172/0,0.14,0.188/0,0.14,0.204/0,0.14,0.22/0,0.14,0.236/0,0.14,0.252/0,0.14,0.268/0,0.14,0.284/1,-0.004,-0.004/1,-0.004,0.012/1,-0.004,0.028/1,-0.004,0.044/1,-0.004,0.06/1,-0.004,0.076/1,-0.004,0.092/1,-0.004,0.108/1,-0.004,0.124/1,-0.004,0.14/1,-0.004,0.156/1,-0.004,0.172/1,-0.004,0.188/1,-0.004,0.204/1,-0.004,0.22/1,-0.004,0.236/1,-0.004,0.252/1,-0.004,0.268/1,-0.004,0.284/1,-0.004,0.3/1,-0.004,0.316/1,0.012,-0.004/1,0.028,-0.004/1,0.044,-0.004/1,0.06,-0.004/1,0.076,-0.004/1,0.092,-0.004/1,0.108,-0.004/1,0.124,-0.004/1,0.14,-0.004/1,0.156,-0.004/1,0.172,-0.004/1,0.188,-0.004/1,0.204,-0.004/1,0.22,-0.004/1,0.236,-0.004/1,0.252,-0.004/1,0.268,-0.004/1,0.284,-0.004/1,0.3,-0.004/1,0.316,-0.004/1,0.332,-0.004/1,0.348,-0.004/1,0.364,-0.004/1,0.38,-0.004/1,0.396,-0.004/1,0.412,-0.004/1,0.428,-0.004/1,0.444,-0.004/1,0.46,-0.004/1,0.476,-0.004/1,0.492,-0.004/1,0.508,-0.004/1,0.524,-0.004/1,0.54,-0.004/1,0.556,-0.004/1,0.572,-0.004/1,0.572,0.012/1,0.572,0.028/1,0.572,0.044/1,0.572,0.06/1,0.572,0.076/1,0.572,0.092/1,0.572,0.108/1,0.572,0.124/1,0.572,0.14/1,0.572,0.156/1,0.572,0.172/1,0.572,0.188/1,0.572,0.204/1,0.572,0.22/1,0.572,0.236/1,0.572,0.252/1,0.572,0.268/1,0.572,0.284/1,0.572,0.3/1,0.572,0.316/2,-0.02,-0.02/2,-0.02,-0.004/2,-0.02,0.012/2,-0.02,0.028/2,-0.02,0.044/2,-0.02,0.06/2,-0.02,0.076/2,-0.02,0.092/2,-0.02,0.108/2,-0.02,0.124/2,-0.02,0.14/2,-0.02,0.156/2,-0.02,0.172/2,-0.02,0.188/2,-0.02,0.204/2,-0.02,0.22/2,-0.02,0.236/2,-0.02,0.252/2,-0.02,0.268/2,-0.02,0.284/2,-0.02,0.3/2,-0.02,0.316/2,-0.004,-0.02/2,0.012,-0.02/2,0.028,-0.02/2,0.044,-0.02/2,0.06,-0.02/2,0.076,-0.02/2,0.092,-0.02/2,0.108,-0.02/2,0.124,-0.02/2,0.14,-0.02/2,0.156,-0.02/2,0.172,-0.02/2,0.188,-0.02/2,0.204,-0.02/2,0.22,-0.02/2,0.236,-0.02/2,0.252,-0.02/2,0.268,-0.02/2,0.284,-0.02/2,0.3,-0.02/2,0.316,-0.02/2,0.332,-0.02/2,0.348,-0.02/2,0.364,-0.02/2,0.38,-0.02/2,0.396,-0.02/2,0.412,-0.02/2,0.428,-0.02/2,0.444,-0.02/2,0.46,-0.02/2,0.476,-0.02/2,0.492,-0.02/2,0.508,-0.02/2,0.524,-0.02/2,0.54,-0.02/2,0.556,-0.02/2,0.572,-0.02/2,0.588,-0.02/2,0.588,-0.004/2,0.588,0.012/2,0.588,0.028/2,0.588,0.044/2,0.588,0.06/2,0.588,0.076/2,0.588,0.092/2,0.588,0.108/2,0.588,0.124/2,0.588,0.14/2,0.588,0.156/2,0.588,0.172/2,0.588,0.188/2,0.588,0.204/2,0.588,0.22/2,0.588,0.236/2,0.588,0.252/2,0.588,0.268/2,0.588,0.284/2,0.588,0.3/2,0.588,0.316";
private var prob02:String="0,0.012,0.012/0,0.012,0.028/0,0.012,0.044/0,0.012,0.06/0,0.012,0.076/0,0.012,0.092/0,0.012,0.108/0,0.012,0.124/0,0.012,0.14/0,0.012,0.156/0,0.012,0.172/0,0.012,0.188/0,0.012,0.204/0,0.012,0.22/0,0.012,0.236/0,0.012,0.252/0,0.012,0.268/0,0.012,0.284/0,0.028,0.012/0,0.028,0.028/0,0.028,0.044/0,0.028,0.06/0,0.028,0.076/0,0.028,0.092/0,0.028,0.108/0,0.028,0.124/0,0.028,0.14/0,0.028,0.156/0,0.028,0.172/0,0.028,0.188/0,0.028,0.204/0,0.028,0.22/0,0.028,0.236/0,0.028,0.252/0,0.028,0.268/0,0.028,0.284/0,0.044,0.012/0,0.044,0.028/0,0.044,0.044/0,0.044,0.06/0,0.044,0.076/0,0.044,0.092/0,0.044,0.108/0,0.044,0.124/0,0.044,0.14/0,0.044,0.156/0,0.044,0.172/0,0.044,0.188/0,0.044,0.204/0,0.044,0.22/0,0.044,0.236/0,0.044,0.252/0,0.044,0.268/0,0.044,0.284/0,0.06,0.012/0,0.06,0.028/0,0.06,0.044/0,0.06,0.06/0,0.06,0.076/0,0.06,0.092/0,0.06,0.108/0,0.06,0.124/0,0.06,0.14/0,0.06,0.156/0,0.06,0.172/0,0.06,0.188/0,0.06,0.204/0,0.06,0.22/0,0.06,0.236/0,0.06,0.252/0,0.06,0.268/0,0.06,0.284/0,0.076,0.012/0,0.076,0.028/0,0.076,0.044/0,0.076,0.06/0,0.076,0.076/0,0.076,0.092/0,0.076,0.108/0,0.076,0.124/0,0.076,0.14/0,0.076,0.156/0,0.076,0.172/0,0.076,0.188/0,0.076,0.204/0,0.076,0.22/0,0.076,0.236/0,0.076,0.252/0,0.076,0.268/0,0.076,0.284/0,0.092,0.012/0,0.092,0.028/0,0.092,0.044/0,0.092,0.06/0,0.092,0.076/0,0.092,0.092/0,0.092,0.108/0,0.092,0.124/0,0.092,0.14/0,0.092,0.156/0,0.092,0.172/0,0.092,0.188/0,0.092,0.204/0,0.092,0.22/0,0.092,0.236/0,0.092,0.252/0,0.092,0.268/0,0.092,0.284/0,0.108,0.012/0,0.108,0.028/0,0.108,0.044/0,0.108,0.06/0,0.108,0.076/0,0.108,0.092/0,0.108,0.108/0,0.108,0.124/0,0.108,0.14/0,0.108,0.156/0,0.108,0.172/0,0.108,0.188/0,0.108,0.204/0,0.108,0.22/0,0.108,0.236/0,0.108,0.252/0,0.108,0.268/0,0.108,0.284/0,0.124,0.012/0,0.124,0.028/0,0.124,0.044/0,0.124,0.06/0,0.124,0.076/0,0.124,0.092/0,0.124,0.108/0,0.124,0.124/0,0.124,0.14/0,0.124,0.156/0,0.124,0.172/0,0.124,0.188/0,0.124,0.204/0,0.124,0.22/0,0.124,0.236/0,0.124,0.252/0,0.124,0.268/0,0.124,0.284/0,0.14,0.012/0,0.14,0.028/0,0.14,0.044/0,0.14,0.06/0,0.14,0.076/0,0.14,0.092/0,0.14,0.108/0,0.14,0.124/0,0.14,0.14/0,0.14,0.156/0,0.14,0.172/0,0.14,0.188/0,0.14,0.204/0,0.14,0.22/0,0.14,0.236/0,0.14,0.252/0,0.14,0.268/0,0.14,0.284/1,-0.004,-0.004/1,-0.004,0.012/1,-0.004,0.028/1,-0.004,0.044/1,-0.004,0.06/1,-0.004,0.076/1,-0.004,0.092/1,-0.004,0.108/1,-0.004,0.124/1,-0.004,0.14/1,-0.004,0.156/1,-0.004,0.172/1,-0.004,0.188/1,-0.004,0.204/1,-0.004,0.22/1,-0.004,0.236/1,-0.004,0.252/1,-0.004,0.268/1,-0.004,0.284/1,-0.004,0.3/1,-0.004,0.316/1,0.012,-0.004/1,0.028,-0.004/1,0.044,-0.004/1,0.06,-0.004/1,0.076,-0.004/1,0.092,-0.004/1,0.108,-0.004/1,0.124,-0.004/1,0.14,-0.004/1,0.156,-0.004/1,0.172,-0.004/1,0.188,-0.004/1,0.204,-0.004/1,0.22,-0.004/1,0.236,-0.004/1,0.252,-0.004/1,0.268,-0.004/1,0.284,-0.004/1,0.3,-0.004/1,0.316,-0.004/1,0.332,-0.004/1,0.348,-0.004/1,0.364,-0.004/1,0.38,-0.004/1,0.396,-0.004/1,0.412,-0.004/1,0.428,-0.004/1,0.444,-0.004/1,0.46,-0.004/1,0.476,-0.004/1,0.492,-0.004/1,0.508,-0.004/1,0.524,-0.004/1,0.54,-0.004/1,0.556,-0.004/1,0.572,-0.004/1,0.572,0.012/1,0.572,0.028/1,0.572,0.044/1,0.572,0.14/1,0.572,0.156/1,0.572,0.172/1,0.572,0.188/1,0.572,0.204/1,0.572,0.22/1,0.572,0.236/1,0.572,0.252/1,0.572,0.268/1,0.572,0.284/1,0.572,0.3/1,0.572,0.316/2,-0.02,-0.02/2,-0.02,-0.004/2,-0.02,0.012/2,-0.02,0.028/2,-0.02,0.044/2,-0.02,0.06/2,-0.02,0.076/2,-0.02,0.092/2,-0.02,0.108/2,-0.02,0.124/2,-0.02,0.14/2,-0.02,0.156/2,-0.02,0.172/2,-0.02,0.188/2,-0.02,0.204/2,-0.02,0.22/2,-0.02,0.236/2,-0.02,0.252/2,-0.02,0.268/2,-0.02,0.284/2,-0.02,0.3/2,-0.02,0.316/2,-0.004,-0.02/2,0.012,-0.02/2,0.028,-0.02/2,0.044,-0.02/2,0.06,-0.02/2,0.076,-0.02/2,0.092,-0.02/2,0.108,-0.02/2,0.124,-0.02/2,0.14,-0.02/2,0.156,-0.02/2,0.172,-0.02/2,0.188,-0.02/2,0.204,-0.02/2,0.22,-0.02/2,0.236,-0.02/2,0.252,-0.02/2,0.268,-0.02/2,0.284,-0.02/2,0.3,-0.02/2,0.316,-0.02/2,0.332,-0.02/2,0.348,-0.02/2,0.364,-0.02/2,0.38,-0.02/2,0.396,-0.02/2,0.412,-0.02/2,0.428,-0.02/2,0.444,-0.02/2,0.46,-0.02/2,0.476,-0.02/2,0.492,-0.02/2,0.508,-0.02/2,0.524,-0.02/2,0.54,-0.02/2,0.556,-0.02/2,0.572,-0.02/2,0.588,-0.02/2,0.588,-0.004/2,0.588,0.012/2,0.588,0.028/2,0.588,0.044/2,0.588,0.14/2,0.588,0.156/2,0.588,0.172/2,0.588,0.188/2,0.588,0.204/2,0.588,0.22/2,0.588,0.236/2,0.588,0.252/2,0.588,0.268/2,0.588,0.284/2,0.588,0.3/2,0.588,0.316";
private var prob03:String="0,0.012,0.012/0,0.012,0.028/0,0.012,0.044/0,0.012,0.06/0,0.012,0.076/0,0.012,0.092/0,0.012,0.108/0,0.012,0.124/0,0.012,0.14/0,0.012,0.156/0,0.012,0.172/0,0.012,0.188/0,0.012,0.204/0,0.012,0.22/0,0.012,0.236/0,0.012,0.252/0,0.012,0.268/0,0.012,0.284/0,0.028,0.012/0,0.028,0.028/0,0.028,0.044/0,0.028,0.06/0,0.028,0.076/0,0.028,0.092/0,0.028,0.108/0,0.028,0.124/0,0.028,0.14/0,0.028,0.156/0,0.028,0.172/0,0.028,0.188/0,0.028,0.204/0,0.028,0.22/0,0.028,0.236/0,0.028,0.252/0,0.028,0.268/0,0.028,0.284/0,0.044,0.012/0,0.044,0.028/0,0.044,0.044/0,0.044,0.06/0,0.044,0.076/0,0.044,0.092/0,0.044,0.108/0,0.044,0.124/0,0.044,0.14/0,0.044,0.156/0,0.044,0.172/0,0.044,0.188/0,0.044,0.204/0,0.044,0.22/0,0.044,0.236/0,0.044,0.252/0,0.044,0.268/0,0.044,0.284/0,0.06,0.012/0,0.06,0.028/0,0.06,0.044/0,0.06,0.06/0,0.06,0.076/0,0.06,0.092/0,0.06,0.108/0,0.06,0.124/0,0.06,0.14/0,0.06,0.156/0,0.06,0.172/0,0.06,0.188/0,0.06,0.204/0,0.06,0.22/0,0.06,0.236/0,0.06,0.252/0,0.06,0.268/0,0.06,0.284/0,0.076,0.012/0,0.076,0.028/0,0.076,0.044/0,0.076,0.06/0,0.076,0.076/0,0.076,0.092/0,0.076,0.108/0,0.076,0.124/0,0.076,0.14/0,0.076,0.156/0,0.076,0.172/0,0.076,0.188/0,0.076,0.204/0,0.076,0.22/0,0.076,0.236/0,0.076,0.252/0,0.076,0.268/0,0.076,0.284/0,0.092,0.012/0,0.092,0.028/0,0.092,0.044/0,0.092,0.06/0,0.092,0.076/0,0.092,0.092/0,0.092,0.108/0,0.092,0.124/0,0.092,0.14/0,0.092,0.156/0,0.092,0.172/0,0.092,0.188/0,0.092,0.204/0,0.092,0.22/0,0.092,0.236/0,0.092,0.252/0,0.092,0.268/0,0.092,0.284/0,0.108,0.012/0,0.108,0.028/0,0.108,0.044/0,0.108,0.06/0,0.108,0.076/0,0.108,0.092/0,0.108,0.108/0,0.108,0.124/0,0.108,0.14/0,0.108,0.156/0,0.108,0.172/0,0.108,0.188/0,0.108,0.204/0,0.108,0.22/0,0.108,0.236/0,0.108,0.252/0,0.108,0.268/0,0.108,0.284/0,0.124,0.012/0,0.124,0.028/0,0.124,0.044/0,0.124,0.06/0,0.124,0.076/0,0.124,0.092/0,0.124,0.108/0,0.124,0.124/0,0.124,0.14/0,0.124,0.156/0,0.124,0.172/0,0.124,0.188/0,0.124,0.204/0,0.124,0.22/0,0.124,0.236/0,0.124,0.252/0,0.124,0.268/0,0.124,0.284/0,0.14,0.012/0,0.14,0.028/0,0.14,0.044/0,0.14,0.06/0,0.14,0.076/0,0.14,0.092/0,0.14,0.108/0,0.14,0.124/0,0.14,0.14/0,0.14,0.156/0,0.14,0.172/0,0.14,0.188/0,0.14,0.204/0,0.14,0.22/0,0.14,0.236/0,0.14,0.252/0,0.14,0.268/0,0.14,0.284/1,-0.004,-0.004/1,-0.004,0.012/1,-0.004,0.028/1,-0.004,0.044/1,-0.004,0.06/1,-0.004,0.076/1,-0.004,0.092/1,-0.004,0.108/1,-0.004,0.124/1,-0.004,0.14/1,-0.004,0.156/1,-0.004,0.172/1,-0.004,0.188/1,-0.004,0.204/1,-0.004,0.22/1,-0.004,0.236/1,-0.004,0.252/1,-0.004,0.268/1,-0.004,0.284/1,-0.004,0.3/1,-0.004,0.316/1,0.012,-0.004/1,0.028,-0.004/1,0.044,-0.004/1,0.06,-0.004/1,0.076,-0.004/1,0.092,-0.004/1,0.108,-0.004/1,0.124,-0.004/1,0.14,-0.004/1,0.156,-0.004/1,0.172,-0.004/1,0.188,-0.004/1,0.204,-0.004/1,0.22,-0.004/1,0.236,-0.004/1,0.252,-0.004/1,0.268,-0.004/1,0.284,-0.004/1,0.3,-0.004/1,0.316,-0.004/1,0.332,-0.004/1,0.348,-0.004/2,0.364,-0.004/1,0.38,-0.004/1,0.396,-0.004/1,0.412,-0.004/1,0.428,-0.004/1,0.444,-0.004/1,0.46,-0.004/1,0.476,-0.004/1,0.492,-0.004/1,0.508,-0.004/1,0.524,-0.004/1,0.54,-0.004/1,0.556,-0.004/1,0.572,-0.004/1,0.572,0.012/1,0.572,0.028/1,0.572,0.044/1,0.572,0.14/1,0.572,0.156/1,0.572,0.172/1,0.572,0.188/1,0.572,0.204/1,0.572,0.22/1,0.572,0.236/1,0.572,0.252/1,0.572,0.268/1,0.572,0.284/1,0.572,0.3/1,0.572,0.316/2,-0.02,-0.02/2,-0.02,-0.004/2,-0.02,0.012/2,-0.02,0.028/2,-0.02,0.044/2,-0.02,0.06/2,-0.02,0.076/2,-0.02,0.092/2,-0.02,0.108/2,-0.02,0.124/2,-0.02,0.14/2,-0.02,0.156/2,-0.02,0.172/2,-0.02,0.188/2,-0.02,0.204/2,-0.02,0.22/2,-0.02,0.236/2,-0.02,0.252/2,-0.02,0.268/2,-0.02,0.284/2,-0.02,0.3/2,-0.02,0.316/2,-0.004,-0.02/2,0.012,-0.02/2,0.028,-0.02/2,0.044,-0.02/2,0.06,-0.02/2,0.076,-0.02/2,0.092,-0.02/2,0.108,-0.02/2,0.124,-0.02/2,0.14,-0.02/2,0.156,-0.02/2,0.172,-0.02/2,0.188,-0.02/2,0.204,-0.02/2,0.22,-0.02/2,0.236,-0.02/2,0.252,-0.02/2,0.268,-0.02/2,0.284,-0.02/2,0.3,-0.02/2,0.316,-0.02/2,0.332,-0.02/2,0.348,-0.02/2,0.364,-0.02/2,0.38,-0.02/2,0.396,-0.02/2,0.412,-0.02/2,0.428,-0.02/2,0.444,-0.02/2,0.46,-0.02/2,0.476,-0.02/2,0.492,-0.02/2,0.508,-0.02/2,0.524,-0.02/2,0.54,-0.02/2,0.556,-0.02/2,0.572,-0.02/2,0.588,-0.02/2,0.588,-0.004/2,0.588,0.012/2,0.588,0.028/2,0.588,0.044/2,0.588,0.14/2,0.588,0.156/2,0.588,0.172/2,0.588,0.188/2,0.588,0.204/2,0.588,0.22/2,0.588,0.236/2,0.588,0.252/2,0.588,0.268/2,0.588,0.284/2,0.588,0.3/2,0.588,0.316/1,0.348,0.012/1,0.348,0.028/1,0.348,0.044/1,0.348,0.06/1,0.348,0.076/1,0.348,0.092/2,0.364,0.012/2,0.364,0.028/2,0.364,0.044/2,0.364,0.06/2,0.364,0.076/1,0.364,0.092/1,0.38,0.012/1,0.38,0.028/1,0.38,0.044/1,0.38,0.06/1,0.38,0.076/1,0.38,0.092";
public function MPS3() {
mc = new MovieClip();
bitmap = new BitmapData(stage.stageWidth, stage.stageHeight, false, 0x000000);
mc.addChild(new Bitmap(bitmap));
addChild(mc);
mc2 = new MovieClip();
addChild(mc2);
cont = new MPSController();
scale = 480.0;
time = new TextField();
var tf:TextFormat = new TextFormat();
tf.color = 0xffffff;
time.defaultTextFormat = tf;
time.text = "time=0.0";
time.x = 10;
time.y = 10;
mc2.addChild(time);
progress=new ProgressBar();
progress.x = 140;
progress.y = 60;
addChild(progress);
rec = new Vector.<Array>();
times = new Vector.<Number>();
var prob1 : RadioButton = new RadioButton(this, 420, 5, "No.1", true, function(e:Event):void { startPb(prob01); } );
var prob2 : RadioButton = new RadioButton(this, 420, 20, "No.2",false, function(e:Event):void { startPb(prob02); } );
var prob3 : RadioButton = new RadioButton(this, 420, 35, "No.3",false, function(e:Event):void { startPb(prob03); } );
addChild(prob1); addChild(prob2); addChild(prob3);
addEventListener(Event.ENTER_FRAME, update);
stage.addEventListener(MouseEvent.MOUSE_DOWN, onDown);
startPb(prob01);
}
private function onDown(e:MouseEvent):void {
counter = 0.0;
isRendering = !isRendering;
}
private function startPb(prob:String):void{
while (cont.array.length > 0) cont.array.pop();
while (rec.length > 0) rec.pop();
while (times.length > 0) times.pop();
isRendering = true
var data:Array = new Array();
var tmp:Array = prob.split("/");
var l:Number = tmp.length;
for (var i:Number = 0; i < l; i++) data.push(tmp[i].split(","));
cont.aveParticleDis = 0.016;
radius = cont.aveParticleDis / 2.0;
ball = new Vector.<Metaball>();
var cs:Vector.<Array> = createPallet();
for (i = 0; i < data.length; i++) {
tmp = data[i];
var p:Particle = new Particle(parseInt(tmp[0]), parseFloat(tmp[1]), parseFloat(tmp[2]), 0, 0, 0, 0);
cont.array.push(p);
var m:Metaball = new Metaball(p, cs, radius, scale);
m.init();
ball.push(m);
}
cont.init();
isThinking = true;
counter = 0;
nextTime = 0.02;
draw();
}
private function update(e:Event):void {
if (isThinking) {
if (cont.totalTime > nextTime) record();
cont.solve();
if (cont.totalTime > endTime) {
record();
isThinking = false;
progress.clear();
}else {
progress.update(cont.totalTime / endTime);
}
}else if (rec.length > counter) {
time.text = "time=" + times[counter].toString();
time.width = time.textWidth + 10;
var i:int = 0, it:int = 0;
var pos:Array = rec[counter];
for each(var p:Particle in cont.array) {
p.x = pos[i++]; p.y = pos[i++];
}
draw();
counter++;
}
}
private function draw():void {
var it:int = 0;
mc2.graphics.clear();
bitmap.fillRect(bitmap.rect, 0x000000);
for each(var p:Particle in cont.array) {
if (p.type == 0 && isRendering) {
ball[it].draw(bitmap);
}else {
mc2.graphics.beginFill(colors[p.type]);
mc2.graphics.drawCircle(p.x * scale + 120, -(p.y * scale-400), radius * scale);
mc2.graphics.endFill();
}
it++;
}
}
private function record():void {
var tmp:Array = new Array();
for each(var p:Particle in cont.array) {
tmp.push(p.x); tmp.push(p.y);
}
rec.push(tmp);
times.push(cont.totalTime);
nextTime += 0.02;
}
private function createPallet():Vector.<Array>{
var ret:Vector.<Array> = new Vector.<Array>();
var r:int, g:int, b:int;
for (var i:int=0;i<256;i++){
r = g = b = 0;
if (i >= 0) b = Math.min(255, 4 * i * 2);
if (i >= 2) g = Math.min(255, 4 * (i / 8));
if (i >= 4) r = Math.min(255, 4 * (i / 16));
ret[i]=[r,g,b];
}
return ret;
}
}
}
import flash.display.MovieClip;
import flash.text.TextField;
import flash.text.TextFormat;
import flash.display.BitmapData;
class ProgressBar extends MovieClip {
private var mes:TextField;
public function ProgressBar() {
super();
mes = new TextField();
var tf:TextFormat = new TextFormat();
tf.color = 0xffffff;
tf.size = 16;
mes.defaultTextFormat = tf;
mes.x = 90;
addChild(mes);
}
public function update(r:Number):void {
graphics.clear();
graphics.lineStyle(1, 0x0000ff);
graphics.drawRect(0, 0, 200, 30);
graphics.beginFill(0x0000ff);
graphics.drawRect(0, 0, 200 * r, 30);
graphics.endFill();
r = Math.round(r * 100.0);
mes.text = r.toString() + "%";
mes.width = mes.textWidth + 5;
}
public function clear():void {
graphics.clear();
mes.text = "";
}
}
class Metaball {
public var colors:Vector.<Array> ;
private var p:Particle;
private var pixels:Vector.<Pixel>;
private var RADIUS:int;
private var NUM_PIXELS:int;
private var radius:Number;
private var scale:Number;
public function Metaball(p:Particle,c:Vector.<Array>,radius:Number,scale:Number){
this.p = p;
this.radius = radius;
this.scale = scale;
RADIUS = Math.floor(radius * scale * 3) as int;
NUM_PIXELS=(2 * RADIUS) * (2 * RADIUS);
pixels = new Vector.<Pixel>();
for (var i:int = 0; i < NUM_PIXELS; i++) pixels[i] = new Pixel();
colors=c;
}
public function init():void {
var no:int;
var c:int = 0;
for (var i:int = -RADIUS; i < RADIUS; i++) {
for (var j:int = -RADIUS; j < RADIUS; j++) {
var z:Number = RADIUS * RADIUS - i * i - j * j;
if (z < 0) {
no = 0;
}else{
z = Math.sqrt(z);
var t:Number = z / RADIUS;
no = (int) (255 * (t * t * t * t));
if (no > 255)no = 255;
if (no < 0)no = 0;
}
pixels[c].dx = i;
pixels[c].dy = j;
pixels[c].no = no;
c++;
}
}
}
public function draw(bitmap:BitmapData):void {
var x:int = Math.round(p.x * scale + 120) as int;
var y:int = Math.round( -(p.y * scale-400)) as int;
for (var i:int = 0; i < NUM_PIXELS; i++) {
var sx:int = x + pixels[i].dx;
if (sx < 0 || sx > bitmap.width-1)continue;
var sy:int = y + pixels[i].dy;
if (sy < 0 || sy > bitmap.height-1)continue;
var no:int = pixels[i].no;
var color:Array = colors[no];
var c:uint = bitmap.getPixel(sx, sy);
var p:Array = getRgb(c);
for (var j:int = 0; j < 3; j++) {
p[j] += color[j];
if (p[j] > 255)p[j] = 255;
}
bitmap.setPixel(sx, sy, get16(p[0], p[1], p[2]));
}
}
private function get16(rr:int,gg:int,bb:int):uint {
return (rr << 16) + (gg << 8) + bb;
}
private function getRgb(c:uint):Array {
var r:int = (c & 0xff0000) >> 16;
var g:int = (c & 0xff00) >> 8;
var b:int = (c & 0xff);
return [r, g, b];
}
}
class Pixel {
public var dx:int;
public var dy:int;
public var no:int;
}
class Particle {
private static const WATER:int = 0;
private static const WALL:int = 2;
private static const PWALL:int = 1;
private static const RHO:Array = [1000.0, 1000.0, 1000.0]; //密度
private static const ALPHA:Array = [0.0000001, 0.0000001, 0.0000001];//比圧縮率
private static var DIM:Number = 2;
private static var neighSize:int = 100 ;
public var press:Number=0.0;
public var dencity:Number=0.0;
public var bcon:int=0 , check:int=0;
public var x:Number=0.0 , y:Number=0.0;
public var vx:Number=0.0 , vy:Number=0.0;
private var dx:Number=0.0 , dy:Number=0.0;
public var poiss:Array , ic:Array;
public var dp:Number=0.0, source:Number=0.0, pmin:Number=0.0;
public var type:int=0 ;
public function Particle( type:int , x:Number , y:Number , vx:Number , vy:Number , p:Number , d:Number ) {
this.x = x ; this.y = y ;
this.type = type ;
press = p ; dencity = d ;
dp = 0.0 ;
this.vx = vx ; this.vy = vy ;
dx = 0.0 ; dy = 0.0 ;
poiss = createArray(neighSize + 1);
ic = createArray(neighSize + 1);
}
private function createArray(n:int):Array {
var ret:Array = new Array();
for (var i:int = 0; i < n; i++) ret[i] = 0;
return ret;
}
public function distanceSq( p:Particle ):Number {
var xx:Number = x - p.x;
var yy:Number = y - p.y;
return ( xx * xx + yy * yy) ;
}
public function getVerocity():Number {
return Math.sqrt(vx * vx + vy * vy) ;
}
public function calcBuoyancy( g:Number , dt:Number ):void {
vy -= g * dt ;
}
public function calcPNDensity( re:Number , array:Array , neigh:Array , num:int ):void {
dencity = calcDensity(re , array , neigh , num) ;
}
public function calcDensity( re:Number , array:Array , neigh:Array , num:int ):Number {
var ret:Number = 0 ;
for( var i:int = 0 ; i < num ; i++ ) {
var r:Number = Math.sqrt(distanceSq(array[neigh[i]]));
var w:Number = getWeight(r , re) ;
ret += w ;
}
return ret ;
}
public static function getWeight( r:Number , re:Number ):Number {
var rr:Number = r / re ;
if( rr < 1.0 ) {
return 1.0 / rr - 1.0 ;
}
else {
return 0.0 ;
}
}
public function calcConvection( dt:Number ):void {
if( isWall() || isPWall() ) return;
x += vx * dt ;
y += vy * dt ;
}
public function calcConvectionDv( dt:Number ):void {
if( isWall() || isPWall() ) return;
x += dx * dt ;
y += dy * dt ;
}
public function isWall():Boolean { return (type == WALL);}
public function isPWall():Boolean { return (type == PWALL);}
public function isWater():Boolean { return (type == WATER);}
public function calcPressureGradientTerm( dt:Number , re:Number , n0p:Number , array:Array , neigh:Array , num:int ):void {
for( var i:int = 0 ; i < num ; i++ ) {
var p:Particle = array[neigh[i]] ;
if( p.isWall()) continue ;
var rr:Number = Math.sqrt(distanceSq(p));
var ddv:Number = dt * ( p.dp - pmin ) / rr ;
ddv = ddv * getWeight(rr , re) / RHO[type];
ddv = ddv / n0p * DIM;
dx += ddv * (x - p.x) / rr ;
dy += ddv * (y - p.y) / rr ;
}
}
public function calcCollision( dt:Number , r2lim:Number , col_rat:Number , id:int , array:Array , neigh:Array , num:int ):void {
var m1:Number = RHO[type] ;
var gx:Number , gy:Number , gz:Number ;
for( var i:int = 0 ; i < num ; i++ ) {
if( neigh[i] <= id ) continue ;
var p:Particle = array[neigh[i]] ;
var xx:Number = p.x - x ;
var yy:Number = p.y - y ;
var rr2:Number = distanceSq(p) ;
if( rr2 < r2lim ) {
var rr:Number = Math.sqrt(rr2) ;
var m2:Number = RHO[p.type] ;
var mm:Number = m1 + m2 ;
gx = ( m1 * vx + m2 * p.vx ) / mm ;
gy = ( m1 * vy + m2 * p.vy ) / mm ;
gx = m1 * ( vx - gx ) ;
gy = m1 * ( vy - gy ) ;
var vabs:Number = ( gx * xx + gy * yy) / rr ;
if( vabs < 0.0 ) continue ;
var vrat:Number = 1.0 + col_rat ;
gx = vrat * vabs * xx / rr ;
gy = vrat * vabs * yy / rr ;
if( !isWall() && !isPWall() ) {
vx -= gx / m1 ;
vy -= gy / m1 ;
x -= dt * gx / m1 ;
y -= dt * gy / m1 ;
}
if( !p.isWall() && !p.isPWall() ) {
p.vx += gx / m2 ;
p.vy += gy / m2 ;
p.x += dt * gx / m2 ;
p.y += dt * gy / m2 ;
}
}
}
}
public function setMatrix( dt:Number , lambda:Number , re:Number , n0:Number , array:Array , iccg:Array , num:int ):void {
poiss[0] = 0 ;
for( var i:int = 0 ; i < num ; i++ ) {
var p:Particle = array[iccg[i]] ;
if( p.bcon == -1 ) {
poiss[i + 1] = 0.0 ;
}else {
var rr:Number = Math.sqrt(distanceSq(p));
var val:Number = 2.0 * DIM / lambda * getWeight(rr , re) / n0 / ((RHO[type]+RHO[p.type])/2.0);
poiss[i + 1] = -val ;
poiss[0] += val ;
}
}
poiss[0] += ALPHA[type] / dt / dt ;
}
public function incom_decomp( id:int , array:Array , iccg:Array , num:Array ):void {
if( bcon != 0 ) return ;
var sum:Number = poiss[0] ;
for( var i:int = 0 ; i < num[id][1] ; i++ ) {
if( iccg[id][i] > id ) continue ;
var p:Particle = array[iccg[id][i]] ;
if( p.bcon != 0 ) continue ;
sum = sum - ic[i + 1] * ic[i + 1] ;
}
ic[0] = Math.sqrt(sum) ;
for(i = 0 ; i < num[id][1] ; i++ ) {
if( iccg[id][i] < id ) continue ;
p = array[iccg[id][i]] ;
if( p.bcon != 0 ) continue ;
sum = poiss[i + 1] ;
var id2:int = iccg[id][i] ;
for( var mj:int = 0 ; mj < num[id2][1] ; mj++ ) {
if( iccg[id2][mj] >= id ) continue ;
var q:Particle = array[iccg[id2][mj]] ;
if( q.bcon != 0 ) continue ;
for( var mi:int = 0 ; mi < num[id][1] ; mi++ ) {
if( iccg[id2][mj] == iccg[id][mi] ) {
sum = sum - ic[mi + 1] * p.ic[mj + 1] ;
break ;
}
}
}
ic[i + 1] = sum / ic[0] ;
for(mj = 0 ; mj < num[id2][1] ; mj++ ) {
if( id == iccg[id2][mj] ) {
p.ic[mj + 1] = ic[i + 1] ;
break ;
}
}
}
}
public function setBcon( nop:Number , cutoff:Number ):void {
if(isWall() ) {
bcon = -1 ; check = -1 ;
}else if( dencity / nop < cutoff ) {
bcon = 1 ; check = 1 ;
}else {
bcon = 0 ; check = 0 ;
}
}
public function checkBcon( array:Array , neigh:Array , num:int ):void {
for( var i:int = 0 ; i < num ; i++ ) {
var p:Particle = array[neigh[i]] ;
if( p.check == 0 ) p.check = 1 ;
}
check = 2 ;
}
public function setSource( dt:Number , nop:Number ):void {
source = 1.0 / dt / dt * ( dencity - nop ) / nop ;
}
public function coeffCheck( diag:Number ):void {
if( check == 0 )poiss[0] *= diag ;
}
public function calcDv( px:Particle , dt:Number , rep:Number , nOp:Number ):void {
var rr:Number = Math.sqrt(distanceSq(px));
var ddv:Number = dt * ( px.dp - pmin ) / rr ;
ddv = ddv * getWeight(rr , rep) / RHO[type];
ddv = ddv / nOp * DIM;
dx += ddv * ( x - px.x ) / rr ;
dy += ddv * ( y - px.y ) / rr ;
}
public function updatePDv():void {
press += dp ;
vx += dx ; vy += dy ;
dp = 0 ;
dx = 0.0 ; dy = 0.0 ;
}
}
class MPSController {
public var aveParticleDis:Number; //初期配置における近接粒子間距離
public const GRAV:Number=9.8; //重力加速度
public const DT_MAX:Number=0.01; //時間増分上限
public const CR_MAX:Number=0.2; //クーラン数上限
public const EPSP:Number=0.000000001; //圧力計算における収束条件
public const IMPX:int=50; //圧力計算における反復数の上限
public const IMNP:int=10; //圧力計算における反復数の下限
public const RE_PND:Number=2.1; //粒子数密度の計算に用いる重み関数の半径
public const RE_ICCG:Number=4.0; //圧力のラプラシアン計算に用いる重み関数の半径
public const MIN_LEN:Number=0.5; //粒子間の再接近距離
public const DIRRICHLET:Number=0.97; //自由表面の条件を与えるパラメータ
public const COL_RATE:Number=0.2; //接近した粒子の反発率
public const LAM_RATE:Number=0.2; //圧力計算の緩和係数
public const INCREASE_DT:Number = 1.1 ;
public var n0p:Number=0.0; //repに対応する粒子数密度
public var n0iccg:Number = 0.0; //reiccgに対応する粒子数密度
public var n0pperticle:int=0; //check
public var array:Array;
public var totalTime:Number ;
private var rep:Number ,rep2:Number, reiccg:Number , reiccg2:Number ;
private var r2lim:Number , ra:Number , lambda:Number ;
private var dt:Number ;
private var DIAG:Number = 2.0;
private var maxNeigh:int = 100;
private var neigh:Array, iccg:Array, num:Array ;
private var r:Array, q:Array;
public function MPSController(){
array = new Array() ;
}
public function init():void {
if ( array.length == 0 ) return ;
num = new Array();
neigh = new Array();
iccg = new Array();
r = new Array();
q = new Array();
for (var i:int = 0; i < array.length; i++) {
num[i] = [0, 0];
var t0:Array = [];
var t1:Array = [];
for (var j:int = 0; j < 100; j++) {
t0[j] = 0.0; t1[j] = 0.0;
}
neigh[i] = t0; iccg[i] = t1;
r[i] = 0.0; q[i] = 0.0;
}
totalTime = 0 ;
rep = aveParticleDis * RE_PND ;
reiccg = aveParticleDis * RE_ICCG ;
reiccg2 = reiccg * reiccg ;
rep2 = rep * rep;
if( n0p == 0.0 ) {
setNeighbor(rep , reiccg) ;
var p:Particle = array[n0pperticle] ;
p.calcPNDensity(rep , array , neigh[n0pperticle] , num[n0pperticle][0]) ;
n0p = p.dencity;
n0iccg = p.calcDensity(reiccg , array , iccg[n0pperticle] , num[n0pperticle][1]) ;
}
r2lim = Math.pow(aveParticleDis , 2) * Math.pow(MIN_LEN , 2) ;
ra = aveParticleDis / 1.7724539 ;
lambda = LAM_RATE * ( reiccg2 * reiccg2 / 12.0 - reiccg * ra * ra * ra / 3.0 + ra * ra * ra * ra / 4.0 ) / ( reiccg2 / 2.0 - reiccg * ra + ra * ra / 2.0 ) ;
dt = DT_MAX;
setNeighbor(rep , reiccg);
for (i = 0 ; i < array.length ; i++) array[i].calcPNDensity(rep , array , neigh[i] , num[i][0]);
}
public function solve():void {
dt = calcDt(CR_MAX , aveParticleDis , DT_MAX , dt) ;
totalTime += dt ;
for( var i:int = 0 ; i < array.length ; i++)array[i].calcBuoyancy(GRAV , dt);
for (i = 0 ; i < array.length ; i++) array[i].calcConvection(dt);
for (i = 0 ; i < array.length ; i++) array[i].calcCollision(dt, r2lim, COL_RATE, i, array, neigh[i], num[i][0]);
setNeighbor(rep , reiccg) ;
for (i = 0 ; i < array.length ; i++ ) array[i].calcPNDensity(rep , array , neigh[i] , num[i][0]);
for(i = 0 ; i < array.length ; i++ )array[i].setBcon(n0p , DIRRICHLET);
setSource(n0p);
checkBcon();
setMatrix();
for(i = 0 ; i < array.length ; i++ )array[i].coeffCheck(DIAG) ;
for(i = 0 ; i < array.length ; i++ )array[i].incom_decomp(i , array , iccg , num) ;
pcg_solver(EPSP, IMNP, IMPX, dt);
for(i = 0 ; i < array.length ; i++ ) if( array[i].dp < 0.0 ) array[i].dp = 0.0 ;
set_pmin() ;
rev_pressgrad(dt , rep , n0p) ;
for(i = 0 ; i < array.length ; i++ )array[i].calcConvectionDv(dt) ;
for(i = 0 ; i < array.length ; i++ )array[i].updatePDv() ;
setNeighbor(rep , reiccg) ;
}
private function setNeighbor( rep:Number , reiccg:Number ):void {
for( var i:int = 0 ; i < array.length ; i++ ) {
num[i][0] = 0; num[i][1] = 0;
}
for(i = 0 ; i < array.length ; i++ ) {
for( var j:int = i + 1 ; j < array.length ; j++ ) {
if( array[i].isWall() && array[j].isWall() ) continue ;
var d:Number = array[i].distanceSq(array[j]) ;
if( d <= rep2 ) {
neigh[i][num[i][0]++] = j ;
neigh[j][num[j][0]++] = i ;
}
if( d <= reiccg2 ) {
iccg[i][num[i][1]++] = j ;
iccg[j][num[j][1]++] = i ;
}
}
}
}
private function setMatrix():void {
for(var i:int = 0 ; i < array.length ; i++ ) {
if( array[i].bcon != 0 ) continue ;
array[i].setMatrix(dt , lambda , reiccg , n0iccg , array , iccg[i] , num[i][1]) ;
}
}
private function calcDt( cr:Number , aveParticleDis:Number , dtmax:Number , dtold:Number ):Number {
var vmax:Number = 0.0 ;
for( var i:int = 0 ; i < array.length ; i++ ) {
if( array[i].isWall() || array[i].isPWall() ) continue ;
var vv:Number = array[i].getVerocity() ;
if( vv > vmax ) vmax = vv ;
}
var dtlimit:Number = dtold * INCREASE_DT ;
if( vmax == 0 ) return dtlimit ;
var dt:Number = cr * aveParticleDis / vmax ;
if( dt > dtlimit ) dt = dtlimit ;
if( dt > dtmax ) dt = dtmax ;
return dt ;
}
private function checkBcon():void {
var count:int,i:int ;
do{
count = 0 ;
for(i = 0 ; i < array.length ; i++ ) {
if( array[i].check == 1 ) {
array[i].checkBcon(array , neigh[i] , num[i][0]) ;
count++ ;
}
}
}while( count != 0 ) ;
}
private function setSource(nop:Number ):void {
var num_dirichlet:int = 0 ;
for( var i:int = 0 ; i < array.length ; i++ ) {
if( array[i].isWall() ) continue ;
if( array[i].bcon == 0 ) {
array[i].setSource(dt , nop) ;
}else if( array[i].bcon == 1 ) {
array[i].source = 0.0;
if( array[i].isWater()) num_dirichlet++ ;
}
}
}
private function pcg_solver( eps:Number , imin:Number , imax:Number , dt:Number ):Boolean {
for( var i:int = 0 ; i < array.length ; i++ ) {
array[i].dp = 0.0 ;
r[i] = 0.0; q[i] = 0.0;
}
var s:Array = mul_matrix_vector() ;
for(i = 0 ; i < array.length ; i++ ) r[i] = array[i].source - s[i] ;
solver_ll(q , r) ;
for(i = 0 ; i < array.length ; i++ )array[i].press = q[i];
var rqo:Number = mul_vector(r , q) ;
for ( var k:int = 0 ; k < imax ; k++ ) {
mul_matrix_vector_press(s) ;
var ps:Number = mul_vector_press(s) ;
var aa:Number = rqo / ps ;
for(i = 0 ; i < array.length ; i++ ) {
var p0:Particle = array[i] ;
p0.dp = p0.dp + aa * p0.press ;
}
for(i = 0 ; i < r.length ; i++ ) r[i] = r[i] - aa * s[i] ;
solver_ll(q , r) ;
var rqn:Number = mul_vector(r , q) ;
var bb:Number = rqn / rqo ;
rqo = rqn ;
for (i = 0 ; i < array.length ; i++ ) array[i].press = (q[i] + bb * array[i].press);
if( checkconv(r , eps , k , imin , imax , dt) ) return true ;
}
return false ;
}
private function mul_matrix_vector():Array {
var y:Array = new Array();
for ( var i:int = 0 ; i < array.length ; i++ ) {
y[i] = 0.0;
if( array[i].bcon != 0 ) continue ;
y[i] = array[i].poiss[0] * array[i].dp ;
for( var j:int = 0 ; j < num[i][1] ; j++ ) {
var p:Particle = array[iccg[i][j]] ;
if( p.bcon != -1 ) continue ;
y[i] = y[i] + p.poiss[j + 1] * p.dp ;
}
}
return y ;
}
private function solver_ll( y:Array , xx:Array ):void {
var x:Array = new Array();
for( var i:int = 0 ; i < xx.length ; i++ ) x[i] = xx[i] ;
for(i = 0 ; i < array.length ; i++ ) {
var p:Particle = array[i] ;
if( p.bcon != 0 ) continue ;
for( var j:int = 0 ; j < num[i][1] ; j++ ) {
var px:Particle = array[iccg[i][j]] ;
if( iccg[i][j] > i || px.bcon != 0 ) continue ;
x[i] = x[i] - p.ic[j + 1] * y[iccg[i][j]] ;
}
y[i] = x[i] / p.ic[0] ;
}
for(i = 0 ; i < y.length ; i++ ) x[i] = y[i] ;
for(i = array.length - 1 ; i >= 0 ; i-- ) {
p = array[i] ;
if( p.bcon != 0 ) continue ;
for(j = 0 ; j < num[i][1] ; j++ ) {
px= array[iccg[i][j]] ;
if( iccg[i][j] < i || px.bcon != 0 ) continue ;
x[i] = x[i] - p.ic[j + 1] * y[iccg[i][j]] ;
}
y[i] = x[i] / p.ic[0] ;
}
}
private function mul_vector( x1:Array , x2:Array ):Number {
var ans:Number = 0.0 ;
for( var i:int = 0 ; i < array.length ; i++ )if( array[i].bcon == 0 ) ans += x1[i] * x2[i] ;
return ans ;
}
private function mul_matrix_vector_press( y:Array ):void {
for ( var i:int = 0 ; i < array.length ; i++ ) {
var p:Particle = array[i] ;
if( p.bcon != 0 ) continue ;
y[i] = p.poiss[0] * p.press ;
for( var j:int = 0 ; j < num[i][1] ; j++ ) {
var px:Particle = array[iccg[i][j]] ;
if( px.bcon == -1 ) continue ;
y[i] = y[i] + p.poiss[j + 1] * px.press;
}
}
}
private function mul_vector_press( x1:Array ):Number {
var ans:Number = 0.0 ;
for ( var i:int = 0 ; i < array.length ; i++ ) if( array[i].bcon == 0 ) ans += x1[i] * array[i].press;
return ans ;
}
private function checkconv( r:Array , eps:Number , k:Number , imin:Number , imax:Number , dt:Number ):Boolean {
var checker:int = 0 ;
var err1sum:Number = 0.0 ;
for( var i:int = 0 ; i < array.length ; i++ ) {
if( array[i].bcon != 0 ) continue ;
var err1:Number = 0.0 ;
var rr:Number = r[i] ;
if ( rr > 0 ) { err1 = rr * dt * dt ; } else { err1 = -rr * dt * dt ; }
if ( err1 > eps ) { checker++; }
err1sum += err1 ;
}
if( checker == 0 && k >= imin ) {
return true ;
}else {
return false ;
}
}
private function set_pmin():void {
for( var i:int = 0 ; i < array.length ; i++ ) {
var p:Particle = array[i] ;
if( p.bcon == -1 || p.isWall() ) continue ;
p.pmin = p.dp ;
for( var j:int = 0 ; j < num[i][1] ; j++ ) {
var px:Particle = array[iccg[i][j]] ;
if( px.isWall() ) continue ;
if( px.dp < p.pmin ) p.pmin = px.dp ;
}
}
}
private function rev_pressgrad( dt:Number , rep:Number , n0p:Number ):void {
for( var i:int = 0 ; i < array.length ; i++ ) {
var p:Particle = array[i] ;
if( p.bcon == -1 || p.isPWall() ) continue ;
for( var j:int = 0 ; j < num[i][1] ; j++ ) {
var px:Particle = array[iccg[i][j]] ;
if( px.isWall() ) continue ;
p.calcDv(px , dt , rep , n0p) ;
}
}
}
}