const FLUID={orificeMaxCm:1.4,atrialFloorCm:2.4,marginCm:{LA_cavity:.3,LV_cavity:.15,RA_cavity:.6,RV_cavity:.3,AORTA:[.35,.2,1.15,.2,.4,.2],PA:[.9,.65,1.05,.65,1.45,.45],SVC:[.45,.2,.25,.2,.85,.3],IVC:[.4,.5,.4,.2,.2,.2],PULM_VEINS:[1.4,.2,1.4,.55,1.5,.2]},ventricles:["LV_cavity","RV_cavity"],link:{LA_cavity:["LV_cavity","MV"],RA_cavity:["RV_cavity","TV"]},vessels:{AORTA:"LV_cavity",PA:"RV_cavity",SVC:"RA_cavity",IVC:"RA_cavity",PULM_VEINS:"LA_cavity"},veinShare:{SVC:.36,IVC:.6,CS:0},freshTau:.6,hops:3,ostiumScale:1.4,capSpeedMax:{valve:340,ostium:150,tip:340},vClamp:350,vortEps:.3},FLUID_PASS_VERT="void main(){ gl_Position = vec4(position.xy, 0.0, 1.0); }";function atlasLayout(A){const[t,s,a]=A,l=Math.ceil(Math.sqrt(a)),i=Math.ceil(a/l);return{nx:t,ny:s,nz:a,cols:l,rows:i,W:l*t,H:i*s}}function fluidMeshVolume(A,t){const s=A.geometry,a=s.attributes.position,l=a.count,i=s.morphAttributes.position||[],e=A.morphTargetDictionary||{},r=new Float32Array(i.length);for(const[v,y]of Object.entries(e))r[y]=t[v]||0;const h=s.morphTargetsRelative;let m=0;for(let v=0;v<r.length;v++)m+=r[v];const n=new Float32Array(l*3);for(let v=0;v<l;v++){let y=a.getX(v),D=a.getY(v),_=a.getZ(v);h||(y*=1-m,D*=1-m,_*=1-m);for(let F=0;F<r.length;F++)r[F]&&(y+=r[F]*i[F].getX(v),D+=r[F]*i[F].getY(v),_+=r[F]*i[F].getZ(v));n[3*v]=y,n[3*v+1]=D,n[3*v+2]=_}const u=s.index?s.index.array:null,o=u?u.length:l;let c=0;for(let v=0;v<o;v+=3){const y=3*(u?u[v]:v),D=3*(u?u[v+1]:v+1),_=3*(u?u[v+2]:v+2);c+=n[y]*(n[D+1]*n[_+2]-n[D+2]*n[_+1])+n[y+1]*(n[D+2]*n[_]-n[D]*n[_+2])+n[y+2]*(n[D]*n[_+1]-n[D+1]*n[_])}return Math.abs(c)/6}class FluidSolver{constructor(t,s,a={}){const l=THREE;this.name=t.name,this.region=t,this.renderer=s,this.h=t.h;const i=a.marginCm||0,e=(Array.isArray(i)?i:[i,i,i,i,i,i]).map(o=>Math.max(0,Math.round(o/this.h)));this.dims=t.dims.map((o,c)=>o+e[c]+e[3+c]),this.origin=t.origin.map((o,c)=>o-e[c]*this.h),this.size=this.dims.map(o=>o*this.h),this.layout=atlasLayout(this.dims);const{W:r,H:h}=this.layout;this.substeps=a.substeps||2,this.cycles=a.cycles||2,this.pre=2,this.post=2,this.coarseIters=a.coarseIters||6,this.maccormack=a.maccormack!==!1,this.vortEps=a.vortEps!=null?a.vortEps:FLUID.vortEps,this.frame=0,this.exportOn=!1,this.dtFrame=0,this.fluxOnly=!1,this.scene=new l.Scene,this.cam=new l.OrthographicCamera(-1,1,1,-1,0,1),this.quad=new l.Mesh(new l.PlaneGeometry(2,2)),this.quad.frustumCulled=!1,this.scene.add(this.quad),this.rts=[];const m=o=>(this.rts.push(o),o);this._rt=(o=r,c=h)=>m(new l.WebGLRenderTarget(o,c,{type:l.FloatType,format:l.RGBAFormat,minFilter:l.NearestFilter,magFilter:l.NearestFilter,depthBuffer:!1,stencilBuffer:!1,generateMipmaps:!1}));const n=o=>{const c=new l.WebGL3DRenderTarget(o[0],o[1],o[2],{type:l.HalfFloatType,format:l.RGBAFormat,minFilter:l.LinearFilter,magFilter:l.LinearFilter,depthBuffer:!1,stencilBuffer:!1,generateMipmaps:!1});return c.texture.wrapS=c.texture.wrapT=c.texture.wrapR=l.ClampToEdgeWrapping,m(c)};for(const o of["tWall","tWallB","tWallVel","tDep","tDepB","tVel","tVelB","tVelC","tVelD","tPhi","tPhiB","tDiv","tVort","tDye","tDyeB","tDyeC","tDyeD"])this[o]=this._rt();this.rW=Math.ceil(r/8),this.rH=Math.ceil(h/8),this.tBlock=this._rt(this.rW,this.rH),this.tBlock2=this._rt(Math.ceil(this.rW/8),Math.ceil(this.rH/8)),this.tStatFrame=this._rt(1,1),this.tStatDiv=this._rt(1,1),this.tCapVel=this._rt(FLUID_MAX_CAPS,1),this.tVol=m(new l.WebGLRenderTarget(r,h,{type:l.HalfFloatType,format:l.RGBAFormat,minFilter:l.LinearFilter,magFilter:l.LinearFilter,depthBuffer:!1,stencilBuffer:!1,generateMipmaps:!1}));const u=t.meta.deform;this.keys=u.keys.filter(o=>t.deform[o]),this.defDims=u.dims,this.w=new Float32Array(7),this.wPrev=new Float32Array(7),this.keys.length>1&&(this.tDispNow=n(u.dims),this.tDispPrev=n(u.dims)),this.dispNowTex=this.dispPrevTex=t.deform[this.keys[0]],this.dispScaleNow=this.dispScalePrev=0,this._initCaps(t),this.tStatCaps=Array.from({length:Math.max(1,Math.ceil(this.caps.length/4))},()=>this._rt(1,1)),this._initMaterials(t),this._initLevels(),this._withState(()=>this.reset())}_initCaps(t){const s=THREE,a=t.insideSection.view,[l,i,e]=t.insideSection.shape,r=t.origin,h=t.h,m=p=>{const f=Math.floor((p.x-r[0])/h),S=Math.floor((p.y-r[1])/h),d=Math.floor((p.z-r[2])/h);return f>=0&&S>=0&&d>=0&&f<l&&S<i&&d<e&&a[f+l*(S+i*d)]>0},n=t.dist.image,[u,o,c]=[n.width,n.height,n.depth],v=p=>{const f=Math.floor((p.x-r[0])/h),S=Math.floor((p.y-r[1])/h),d=Math.floor((p.z-r[2])/h);return f>=0&&S>=0&&d>=0&&f<u&&S<o&&d<c&&n.data[f+u*(S+o*d)]/255*2.5-.5>0},y=(p,f,S)=>p.clone().addScaledVector(f,S),D=(p,f)=>{let S=0;for(let d=1;d<=40;d++)m(y(p,f,d*.075))&&S++;return S},_=(p,f)=>{for(let S=60;S>=-60;S--){const d=S*.075;if(v(y(p,f,d))&&v(y(p,f,d-.075)))return d}return 0};this.caps=[];for(const p of["inlets","outlets"])for(const f of t.meta.bounds[p]||[]){if(this.caps.length>=FLUID_MAX_CAPS||!["valve","ostium","tip"].includes(f.kind))continue;const S=!!FLUID.link[this.name],d={kind:f.kind,valve:f.valve||null,vessel:f.vessel||null,side:p,c:new s.Vector3(...f.center),n:new s.Vector3(...f.axis).normalize(),r:f.kind==="valve"?Math.min(f.radius,FLUID.orificeMaxCm):f.radius*FLUID.ostiumScale};d.rFwd=f.kind==="valve"&&S?Math.min(f.radius,FLUID.atrialFloorCm):d.r,f.kind==="valve"?d.out=d.n.clone().multiplyScalar(p==="inlets"?-1:1):d.out=D(d.c,d.n)<=D(d.c,d.n.clone().negate())?d.n.clone():d.n.clone().negate(),d.uEnd=_(d.c,d.out),this.caps.push(d)}this.capOpen=new Float32Array(FLUID_MAX_CAPS),this.capR=new Float32Array(FLUID_MAX_CAPS),this.caps.forEach((p,f)=>{this.capOpen[f]=p.kind==="valve"?0:1,this.capR[f]=p.rFwd}),this.alpha=new Float32Array(FLUID_MAX_CAPS),this.beta=new Float32Array(FLUID_MAX_CAPS),this.capVmax=new Float32Array(FLUID_MAX_CAPS),this.gamma=new Float32Array(FLUID_MAX_CAPS),this.linkIdx=new Float32Array(FLUID_MAX_CAPS),this.group=new Float32Array(FLUID_MAX_CAPS).fill(-1),this.localLink=new Float32Array(FLUID_MAX_CAPS).fill(-1);const F=()=>Array.from({length:FLUID_MAX_CAPS},()=>new s.Vector4);this.discU=F(),this.axisU=F(),this.endU=F(),this.caps.forEach((p,f)=>{this.capVmax[f]=FLUID.capSpeedMax[p.kind],this.discU[f].set(p.c.x,p.c.y,p.c.z,p.r),this.axisU[f].set(p.n.x,p.n.y,p.n.z,0),this.endU[f].set(p.out.x,p.out.y,p.out.z,p.uEnd)})}_initMaterials(t){const s=THREE,a=t.meta.deform,l=()=>({uDist:{value:t.dist},uDistOrigin:{value:new s.Vector3(...t.origin)},uDistSize:{value:new s.Vector3(...t.size)},uOrigin:{value:new s.Vector3(...this.origin)},uSize:{value:new s.Vector3(...this.size)},uDispNow:{value:null},uDispPrev:{value:null},uDispScaleNow:{value:1},uDispScalePrev:{value:1},uDefOrigin:{value:new s.Vector3(...a.origin)},uDefSize:{value:new s.Vector3(...a.dims.map(n=>n*a.h))},uDims:{value:new s.Vector3(...this.dims)},uCols:{value:this.layout.cols},uH:{value:this.h},uDt:{value:0},tVel:{value:null},tPhi:{value:null},tDiv:{value:null},tWall:{value:null},tWallVel:{value:null},tDep:{value:null},tCapVel:{value:null}}),i=(n,u)=>new s.ShaderMaterial({uniforms:u,vertexShader:FLUID_PASS_VERT,fragmentShader:n,depthTest:!1,depthWrite:!1}),e=(n,u={})=>i(n,{...l(),...u}),r={value:this.caps.length},h={};for(let n=0;n<7;n++)h["uDef"+n]={value:t.deform[this.keys[n]]||t.deform[this.keys[0]]};this.matDisp=i(FLUID_DISP_GLSL,{...h,uW:{value:new Float32Array(7)},uNDef:{value:this.keys.length},uLayer:{value:0}}),this.matWall=e(FLUID_WALL_GLSL,{uDisc:{value:this.discU},uDiscAxis:{value:this.axisU},uDiscEnd:{value:this.endU},uNDisc:r,uFirst:{value:1}}),this.matWallFix=e(FLUID_WALLFIX_GLSL,{uDiscEnd:{value:this.endU}}),this.matWallVel=e(FLUID_WALLVEL_GLSL),this.matDeposit=e(FLUID_DEPOSIT_GLSL,{uHop:{value:0}}),this.matRFrame=e(FLUID_RFRAME_GLSL),this.matRCaps=e(FLUID_RCAPS_GLSL,{uCapBase:{value:0},uDiscEnd:{value:this.endU}}),this.matRDiv=e(FLUID_RDIV_GLSL),this.matSum8=i(FLUID_SUM8_GLSL,{tBlock:{value:null}}),this.matFinal=i(FLUID_FINAL_GLSL,{tBlock:{value:null},uBlockSize:{value:new s.Vector2(this.tBlock2.width,this.tBlock2.height)}}),this.matCapVel=i(FLUID_CAPVEL_GLSL,{tStatFrame:{value:null},tStatCaps0:{value:null},tStatCaps1:{value:null},tStatCaps2:{value:null},tStatCaps3:{value:null},tStatLink:{value:null},tLinkCap:{value:null},uAlpha:{value:this.alpha},uBeta:{value:this.beta},uGamma:{value:this.gamma},uLinkIdx:{value:this.linkIdx},uGroup:{value:this.group},uLocalLink:{value:this.localLink},uCapVmax:{value:this.capVmax},uDiscEnd:{value:this.endU},uNDisc:r,uH:{value:this.h},uDtFrame:{value:0}}),this.matAdvect=e(FLUID_ADVECT_GLSL,{tAdvSrc:{value:null},uAdvDt:{value:0},uVClamp:{value:FLUID.vClamp}}),this.matMac=e(FLUID_MACCORMACK_GLSL,{tAdvSrc:{value:null},tMcBack:{value:null},uAdvDt:{value:0},uVClamp:{value:FLUID.vClamp}}),this.matVort=e(FLUID_VORT_GLSL),this.matConfine=e(FLUID_CONFINE_GLSL,{tVort:{value:null},uVortEps:{value:0}}),this.matDiv=e(FLUID_DIV_GLSL,{uDtFrame:{value:0}}),this.matProject=e(FLUID_PROJECT_GLSL),this.matDyeAdv=e(FLUID_DYE_ADV_GLSL,{tDyeSrc:{value:null},uAdvDt:{value:0}}),this.matDyeFix=e(FLUID_DYE_FIX_GLSL,{tDyeN:{value:null},tDyeHat:{value:null},tDyeBack:{value:null},uAdvDt:{value:0},uDecay:{value:1},uMac:{value:1}}),this.matExport=e(FLUID_EXPORT_GLSL,{tDyeCur:{value:null},tVortCur:{value:null}}),this.matFill=i(FLUID_FILL_GLSL,{uFill:{value:new s.Vector4}});const m=()=>({tTypeF:{value:null},tRhsF:{value:null},tPhiF:{value:null},tResF:{value:null},tTypeC:{value:null},tPhiC:{value:null},tStatDiv:{value:null},uDimsF:{value:new s.Vector3},uDimsC:{value:new s.Vector3},uColsF:{value:1},uColsC:{value:1},uParity:{value:0},uH:{value:this.h},uOmega:{value:1},uUseStatDiv:{value:0}});this.mg={};for(const[n,u]of Object.entries({smooth:FLUID_MG_SMOOTH_GLSL,residual:FLUID_MG_RESIDUAL_GLSL,coarsen:FLUID_MG_COARSEN_GLSL,restrict:FLUID_MG_RESTRICT_GLSL,prolong:FLUID_MG_PROLONG_GLSL}))this.mg[n]=i(u,m());this.mats=[this.matDisp,this.matWall,this.matWallFix,this.matWallVel,this.matDeposit,this.matRFrame,this.matRCaps,this.matRDiv,this.matSum8,this.matFinal,this.matCapVel,this.matAdvect,this.matMac,this.matVort,this.matConfine,this.matDiv,this.matProject,this.matDyeAdv,this.matDyeFix,this.matExport,this.matFill,...Object.values(this.mg)]}_initLevels(){this.levels=[{dims:this.dims,layout:this.layout,h:this.h,tRes:this._rt()}];for(let t=1;t<3;t++){const s=this.levels[t-1],a=s.dims.map(e=>Math.ceil(e/2));if(Math.min(...a)<3)break;const l=atlasLayout(a),i=()=>this._rt(l.W,l.H);this.levels.push({dims:a,layout:l,h:s.h*2,tType:i(),tRhs:i(),tPhi:i(),tPhiB:i(),tRes:i()})}}_withState(t){const s=this.renderer,a={target:s.getRenderTarget(),autoClear:s.autoClear,shadow:s.shadowMap.autoUpdate,xr:s.xr.enabled};s.autoClear=!1,s.shadowMap.autoUpdate=!1,s.xr.enabled=!1;try{t()}finally{s.setRenderTarget(a.target),s.autoClear=a.autoClear,s.shadowMap.autoUpdate=a.shadow,s.xr.enabled=a.xr}}_draw(t,s,a=0){this.quad.material=t,this.renderer.setRenderTarget(s,a),this.renderer.render(this.scene,this.cam)}_run(t,s,a=0,l=0){const i=t.uniforms;i.tWall&&(i.uDt.value=a,i.tVel.value=this.tVel.texture,i.tPhi.value=this.tPhi.texture,i.tDiv.value=this.tDiv.texture,i.tWall.value=this.tWall.texture,i.tWallVel.value=this.tWallVel.texture,i.tDep.value=this.tDep.texture,i.tCapVel.value=this.tCapVel.texture,i.uDispNow.value=this.dispNowTex,i.uDispPrev.value=this.dispPrevTex,i.uDispScaleNow.value=this.dispScaleNow,i.uDispScalePrev.value=this.dispScalePrev),this._draw(t,s,l)}_swap(t,s){const a=this[t];this[t]=this[s],this[s]=a}_fill(t,s){this.matFill.uniforms.uFill.value.set(s[0],s[1],s[2],s[3]),this._draw(this.matFill,t)}_reduce(t,s){this._run(t,this.tBlock),this.matSum8.uniforms.tBlock.value=this.tBlock.texture,this._draw(this.matSum8,this.tBlock2),this.matFinal.uniforms.tBlock.value=this.tBlock2.texture,this._draw(this.matFinal,s)}_disp(t,s){const a=this.matDisp.uniforms;a.uW.value.set(s);for(let l=0;l<this.defDims[2];l++)a.uLayer.value=l,this._draw(this.matDisp,t,l)}reset(){this.resetFlow();for(const t of["tWallB","tDep","tDepB","tCapVel","tStatDiv","tStatFrame"])this._fill(this[t],[0,0,0,0]);for(const t of this.tStatCaps)this._fill(t,[0,0,0,0]);this.frame=0}resetFlow(){for(const t of["tVel","tVelB","tVelC","tVelD","tVort","tDye","tDyeB","tDyeC","tDyeD","tPhi","tPhiB","tDiv"])this._fill(this[t],[0,0,0,0]);for(let t=1;t<this.levels.length;t++)for(const s of["tPhi","tPhiB","tRhs"])this._fill(this.levels[t][s],[0,0,0,0])}setFluxOnly(t){this.fluxOnly&&!t&&this._withState(()=>this.resetFlow()),this.fluxOnly=!!t}step(t,s,a=null){this._withState(()=>this._step(t,s,a))}_step(t,s,a){const l=this.frame===0;for(let e=0;e<7;e++)this.wPrev[e]=this.w[e],this.w[e]=e<this.keys.length&&s[this.keys[e]]||0;l&&this.wPrev.set(this.w),this.tDispNow?(this._swap("tDispNow","tDispPrev"),this._disp(this.tDispNow,this.w),l&&this._disp(this.tDispPrev,this.w),this.dispNowTex=this.tDispNow.texture,this.dispPrevTex=this.tDispPrev.texture,this.dispScaleNow=1,this.dispScalePrev=1):(this.dispNowTex=this.dispPrevTex=this.region.deform[this.keys[0]],this.dispScaleNow=this.w[0],this.dispScalePrev=this.wPrev[0]);for(let e=0;e<this.caps.length;e++)this.axisU[e].w=this.capOpen[e],this.discU[e].w=this.capR[e];this.matWall.uniforms.uFirst.value=l?1:0,this._run(this.matWall,this.tWallB,t),this._swap("tWall","tWallB");for(let e=0;e<2;e++)this._run(this.matWallFix,this.tWallB,t),this._swap("tWall","tWallB");this._run(this.matWallVel,this.tWallVel,t);for(let e=0;e<FLUID.hops;e++)this.matDeposit.uniforms.uHop.value=e,this._run(this.matDeposit,this.tDepB,t),this._swap("tDep","tDepB");this._reduce(this.matRFrame,this.tStatFrame);for(let e=0;e<this.tStatCaps.length;e++)this.matRCaps.uniforms.uCapBase.value=4*e,this._reduce(this.matRCaps,this.tStatCaps[e]);const i=this.matCapVel.uniforms;i.tStatFrame.value=this.tStatFrame.texture;for(let e=0;e<4;e++)i["tStatCaps"+e].value=this.tStatCaps[Math.min(e,this.tStatCaps.length-1)].texture;if(i.tLinkCap.value=a?a.tCapVel.texture:this.tStatFrame.texture,i.tStatLink.value=(a||this).tStatFrame.texture,i.uDtFrame.value=t,this._draw(this.matCapVel,this.tCapVel),this.dtFrame=t,this.frame++,!this.fluxOnly){if(t>0){const e=t/this.substeps;for(let r=0;r<this.substeps;r++)this._substep(e,t);this._dye(t)}this.exportOn&&this._export()}}_substep(t,s){const a=this.matAdvect.uniforms;if(a.tAdvSrc.value=this.tVel.texture,a.uAdvDt.value=t,this._run(this.matAdvect,this.tVelB,t),this.maccormack){a.tAdvSrc.value=this.tVelB.texture,a.uAdvDt.value=-t,this._run(this.matAdvect,this.tVelC,t);const l=this.matMac.uniforms;l.tAdvSrc.value=this.tVelB.texture,l.tMcBack.value=this.tVelC.texture,l.uAdvDt.value=t,this._run(this.matMac,this.tVelD,t),this._swap("tVel","tVelD")}else this._swap("tVel","tVelB");this.vortEps>0&&(this._run(this.matVort,this.tVort,t),this.matConfine.uniforms.tVort.value=this.tVort.texture,this.matConfine.uniforms.uVortEps.value=this.vortEps,this._run(this.matConfine,this.tVelB,t),this._swap("tVel","tVelB")),this.matDiv.uniforms.uDtFrame.value=s,this._run(this.matDiv,this.tDiv,t),this._reduce(this.matRDiv,this.tStatDiv),this._solvePressure(),this._run(this.matProject,this.tVelB,t),this._swap("tVel","tVelB")}_dye(t){const s=this.matDyeAdv.uniforms;s.tDyeSrc.value=this.tDye.texture,s.uAdvDt.value=t,this._run(this.matDyeAdv,this.tDyeB,t),this.maccormack&&(s.tDyeSrc.value=this.tDyeB.texture,s.uAdvDt.value=-t,this._run(this.matDyeAdv,this.tDyeC,t));const a=this.matDyeFix.uniforms;a.tDyeN.value=this.tDye.texture,a.tDyeHat.value=this.tDyeB.texture,a.tDyeBack.value=this.tDyeC.texture,a.uAdvDt.value=t,a.uDecay.value=Math.exp(-t/FLUID.freshTau),a.uMac.value=this.maccormack?1:0,this._run(this.matDyeFix,this.tDyeD,t),this._swap("tDye","tDyeD")}_export(){const t=this.matExport.uniforms;t.tDyeCur.value=this.tDye.texture,t.tVortCur.value=this.tVort.texture,this._run(this.matExport,this.tVol)}_capComponents(){const t=this.region,s=t.insideSection.view,[a,l,i]=t.insideSection.shape,e=t.origin,r=t.h,h=a*l*i,m=new Int32Array(h).fill(-1),n=new Int32Array(h),u=[];for(let c=0;c<h;c++){if(s[c]<=0||m[c]>=0)continue;const v=u.length;let y=0,D=0;for(n[y++]=c,m[c]=v;y;){const _=n[--y];D++;const F=_%a,p=Math.floor(_/a)%l,f=Math.floor(_/(a*l)),S=[F>0?_-1:-1,F<a-1?_+1:-1,p>0?_-a:-1,p<l-1?_+a:-1,f>0?_-a*l:-1,f<i-1?_+a*l:-1];for(const d of S)d>=0&&s[d]>0&&m[d]<0&&(m[d]=v,n[y++]=d)}u.push(D)}const o=c=>{const v=Math.floor((c.x-e[0])/r),y=Math.floor((c.y-e[1])/r),D=Math.floor((c.z-e[2])/r);return v>=0&&y>=0&&D>=0&&v<a&&y<l&&D<i?m[v+a*(y+l*D)]:-1};this.capComp=this.caps.map(c=>{for(let v=.15;v<=2;v+=.15){const y=o(c.c.clone().addScaledVector(c.out,-v));if(y>=0)return y}return-1}),this.compSize=u}bindVol(t){t.uVol.value=this.tVol.texture,t.uVolDims.value.set(...this.dims),t.uVolAtlas.value.set(this.layout.W,this.layout.H),t.uVolCols.value=this.layout.cols}_lv(t){return t===0?{...this.levels[0],tType:this.tWall,tRhs:this.tDiv,tPhi:this.tPhi,tPhiB:this.tPhiB}:this.levels[t]}_swapPhi(t){if(t===0)this._swap("tPhi","tPhiB");else{const s=this.levels[t];[s.tPhi,s.tPhiB]=[s.tPhiB,s.tPhi]}}_mgRun(t,s,a,l=1,i=0){const e=this.mg[t],r=e.uniforms,h=this._lv(a),m=a+1<this.levels.length?this._lv(a+1):null;r.tTypeF.value=h.tType.texture,r.tRhsF.value=h.tRhs.texture,r.tPhiF.value=h.tPhi.texture,r.tResF.value=h.tRes.texture,r.uDimsF.value.set(...h.dims),r.uColsF.value=h.layout.cols,r.uH.value=h.h,m&&(r.tTypeC.value=m.tType.texture,r.tPhiC.value=m.tPhi.texture,r.uDimsC.value.set(...m.dims),r.uColsC.value=m.layout.cols),r.tStatDiv.value=this.tStatDiv.texture,r.uUseStatDiv.value=a===0?1:0,r.uOmega.value=l,r.uParity.value=i,this._draw(e,s)}_smooth(t,s,a){for(let l=0;l<s;l++)for(const i of[0,1])this._mgRun("smooth",this._lv(t).tPhiB,t,a,i),this._swapPhi(t)}_vcycle(t){if(t===this.levels.length-1){this._smooth(t,this.coarseIters,1.5);return}this._smooth(t,this.pre,1),this._mgRun("residual",this._lv(t).tRes,t),this._mgRun("restrict",this._lv(t+1).tRhs,t),this._fill(this._lv(t+1).tPhi,[0,0,0,0]),this._vcycle(t+1),this._mgRun("prolong",this._lv(t).tPhiB,t),this._swapPhi(t),this._smooth(t,this.post,1)}_solvePressure(){for(let t=0;t+1<this.levels.length;t++)this._mgRun("coarsen",this._lv(t+1).tType,t);for(let t=0;t<this.cycles;t++)this._vcycle(0)}_read(t,s=this.layout.W,a=this.layout.H){const l=new Float32Array(s*a*4);return this.renderer.readRenderTargetPixels(t,0,0,s,a,l),l}readStatTex(t){return Array.from(this._read(this[t],t==="tCapVel"?FLUID_MAX_CAPS:1,1))}readStats(){this._withState(()=>{this.matDiv.uniforms.uDtFrame.value=this.dtFrame,this._run(this.matDiv,this.tDiv,0)});const{nx:t,ny:s,nz:a,cols:l,W:i}=this.layout,e=this._read(this.tWall),r=this._read(this.tVel),h=this._read(this.tDiv),m=this._read(this.tDep),n=this._read(this.tDye),u=this._read(this.tCapVel,FLUID_MAX_CAPS,1),o=(w,L,C)=>4*((Math.floor(C/l)*s+L)*i+C%l*t+w),c=this.h,v=c*c,y=v*c;let D=0,_=0,F=0,p=0,f=0,S=0,d=0,b=0,G=0,P=0;const M=new Float64Array(FLUID_MAX_CAPS),N=new Int32Array(FLUID_MAX_CAPS),g=new Float64Array(FLUID_MAX_CAPS),z=w=>w>3.5&&w<19.5;for(let w=0;w<a;w++)for(let L=0;L<s;L++)for(let C=0;C<t;C++){const x=o(C,L,w),I=e[x];if(z(I)){N[Math.round(I-4)]++,d+=e[x+2];continue}if(I>.5)continue;const R=r[x],k=r[x+1],U=r[x+2],V=h[x];if(!Number.isFinite(R+k+U+V)){S++;continue}D++,_+=V*V,F+=V,d+=e[x+2],b+=e[x+2]-e[x+3]-m[x],G+=n[x],P=Math.max(P,n[x]);const B=Math.hypot(R,k,U);p=Math.max(p,B),f+=B;const q=[[1,0,0,()=>-r[o(C+1,L,w)]],[-1,0,0,()=>R],[0,1,0,()=>-r[o(C,L+1,w)+1]],[0,-1,0,()=>k],[0,0,1,()=>-r[o(C,L,w+1)+2]],[0,0,-1,()=>U]];for(const[$,Y,Z,J]of q){const E=C+$,O=L+Y,T=w+Z;if(E<0||O<0||T<0||E>=t||O>=s||T>=a)continue;const j=e[o(E,O,T)];if(!z(j))continue;const W=Math.round(j-4);g[W]=Math.max(g[W],B),e[x+2]>=.5&&(M[W]+=J()*v)}}const X=this.caps.map((w,L)=>({id:w.valve||w.vessel,open:+this.capOpen[L].toFixed(2),cells:N[L],fluxIn:+M[L].toFixed(1),phi:+u[4*L+3].toFixed(1),v:+Math.hypot(u[4*L],u[4*L+1],u[4*L+2]).toFixed(1),vmaxNear:+g[L].toFixed(0)})),H=this.dtFrame||1;return{fluid:D,volumeML:+(d*y).toFixed(2),divRms:D?Math.sqrt(_/D):0,divMean:D?F/D:0,vmax:p,vmean:D?f/D:0,nan:S,qWall:+(b*y/H).toFixed(1),capFluxSum:+M.reduce((w,L)=>w+L,0).toFixed(1),dyeMean:D?G/D:0,dyeMax:P,caps:X}}dispose(){for(const t of this.rts)t.dispose();for(const t of this.mats)t.dispose();this.quad.geometry.dispose()}}const FLUID_REGIONS=["LV_cavity","RV_cavity","LA_cavity","RA_cavity","AORTA","PA","SVC","IVC","PULM_VEINS"];class FluidSystem{constructor(t,s,a={}){this.renderer=s,this.solvers={};const l=(a.chambers||FLUID_REGIONS).filter(e=>t.regions[e]);for(const e of l)this.solvers[e]=new FluidSolver(t.regions[e],s,{...a,marginCm:FLUID.marginCm[e]||0});const i=e=>FLUID.ventricles.includes(e)?0:FLUID.link[e]?1:2;this.order=l.slice().sort((e,r)=>i(e)-i(r)),this.shares={},this.enabled=!0;for(const e of l)FLUID.link[e]&&this._veinFractions(this.solvers[e]);for(const e of l)FLUID.vessels[e]&&this._linkVessel(this.solvers[e])}_linkVessel(t){const s=this.solvers[FLUID.vessels[t.name]]||null;if(t.linkSolver=s,t.linkOf=t.caps.map(()=>-1),!s)return;t.caps.forEach((i,e)=>{if(i.kind==="tip")return;let r=-1,h=1e9;s.caps.forEach((m,n)=>{if(m.kind!==i.kind||(m.valve||m.vessel)!==(i.valve||i.vessel))return;const u=m.c.distanceTo(i.c);u<h&&(h=u,r=n)}),r>=0&&h<1&&(t.linkOf[e]=r)}),t._capComponents();const l=[...new Set(t.capComp.filter(i=>i>=0))].reduce((i,e)=>i+t.compSize[e],0)||1;t.tipPlan=t.caps.map((i,e)=>{if(i.kind!=="tip")return null;const r=t.capComp[e],h=t.caps.map((m,n)=>n).filter(m=>t.linkOf[m]>=0&&t.capComp[m]===r);return{volFrac:r>=0?t.compSize[r]/l:0,local:h.length===1?h[0]:-1,nLinks:h.length}}),console.info("[fluid] "+t.name+" \u63A5 "+s.name+"\uFF1A",t.caps.map((i,e)=>`${i.kind}:${i.valve||i.vessel||"\u65AD\u53E3"}\u2192\u5757${t.capComp[e]}`+(t.linkOf[e]>=0?`=\u5FC3\u8154\u53E3${t.linkOf[e]}`:t.tipPlan[e]?`(${t.tipPlan[e].nLinks} \u4E2A\u63A5\u53E3)`:"")).join(" "))}_veinFractions(t){t.veinFrac=new Float32Array(FLUID_MAX_CAPS);const s={};t.caps.forEach((i,e)=>{i.kind==="ostium"&&(s[i.vessel]=s[i.vessel]||[]).push(e)});const a=Object.fromEntries(Object.keys(s).map(i=>[i,FLUID.veinShare[i]!=null?FLUID.veinShare[i]:1])),l=Object.values(a).reduce((i,e)=>i+e,0)||1;Object.entries(s).forEach(([i,e],r)=>{for(const h of e)t.veinFrac[h]=a[i]/l,t.group[h]=r})}_ventricleShares(t,s,a,l){const i=new Float32Array(FLUID_MAX_CAPS);if(!a||!l)return i;const e=[];if(t.caps.forEach((n,u)=>{if(n.kind!=="valve"||t.capR[u]*Math.sqrt(t.capOpen[u])<1.5*t.h)return;const o=flowAt(a,n.valve,s);e.push([u,n.side==="inlets"?o:-o])}),!e.length)return i;const r=e.reduce((n,[,u])=>n+Math.max(0,l*u),0);if(r>1e-6){for(const[n,u]of e)i[n]=Math.max(0,l*u)/r;return i}const h=e.filter(([n])=>t.caps[n].side==="inlets"==l>0).map(([n])=>[n,t.capR[n]**2*t.capOpen[n]]),m=h.reduce((n,[,u])=>n+u,0);if(m>1e-9)for(const[n,u]of h)i[n]=u/m;return i}step(t,s,a,l,i,e=null){if(this.enabled)for(const r of this.order){const h=this.solvers[r],m=i[r]||{};h.caps.forEach((u,o)=>{h.capOpen[o]=u.kind==="valve"?Math.max(0,Math.min(1,e?e(u.valve):valveOpenAt(t,a,u.valve))):1}),h.alpha.fill(0),h.beta.fill(0),h.gamma.fill(0),h.linkIdx.fill(0);let n=null;if(FLUID.vessels[r])n=h.linkSolver,n?h.caps.forEach((u,o)=>{if(h.linkOf[o]>=0)h.gamma[o]=-1,h.linkIdx[o]=h.linkOf[o];else if(h.tipPlan[o]){const c=h.tipPlan[o];h.alpha[o]=c.volFrac,h.group[o]=10+h.capComp[o],h.localLink[o]=c.local}}):h.caps.forEach((u,o)=>{h.tipPlan&&h.tipPlan[o]&&(h.alpha[o]=h.tipPlan[o].volFrac,h.group[o]=10+h.capComp[o])});else if(FLUID.ventricles.includes(r)){const u=m.V||0,o=h._wvPrev!=null?u-h._wvPrev:0;h._wvPrev=u;const c=Math.abs(o)>1e-7?-Math.sign(o):0;this.shares[r]=this._ventricleShares(h,t,l,c),h.alpha.set(this.shares[r])}else{const[u,o]=FLUID.link[r]||[],c=this.solvers[u],v=c?c.caps.findIndex(_=>_.valve===o):-1,y=c&&v>=0&&this.shares[u]?this.shares[u][v]:0;n=c||null;const D=l?flowAt(l,o,t)<0:!1;h.caps.forEach((_,F)=>{_.kind==="valve"&&(h.capR[F]=D?_.r:_.rFwd),_.kind==="valve"&&_.valve===o?h.beta[F]=c?-y:0:_.kind==="ostium"&&(h.alpha[F]=h.veinFrac[F],h.beta[F]=c?h.veinFrac[F]*y:0)})}h.step(s,m,n)}}setExport(t){for(const s of Object.values(this.solvers))s.exportOn=!!t}setFluxOnly(t){for(const s of Object.values(this.solvers))s.setFluxOnly(t(s.name))}reset(){for(const t of Object.values(this.solvers))t._withState(()=>t.reset()),t._wvPrev=null}dispose(){for(const t of Object.values(this.solvers))t.dispose()}}Object.assign(window,{FluidSolver,FluidSystem,atlasLayout,fluidMeshVolume,FLUID,FLUID_REGIONS});
