{const r=this.nop[n][i]-1;o(` - Applied constant temperature to node ${r+1} (element ${n+1}, local node ${i+1})`),e[r]=s;for(let n=0;n{if("constantTemp"===this.boundaryConditions[n][0]){const s=this.boundaryConditions[n][1];o(`Boundary ${n}: Applying constant temperature of ${s} K (Dirichlet condition)`),this.boundaryElements[n].forEach((([n,i])=>{if("linear"===this.elementOrder){({0:[0],1:[1]})[i].forEach((i=>{const r=this.nop[n][i]-1;o(` - Applied constant temperature to node ${r+1} (element ${n+1}, local node ${i+1})`),e[r]=1,t[r]=s}))}else if("quadratic"===this.elementOrder){({0:[0],2:[2]})[i].forEach((i=>{const r=this.nop[n][i]-1;o(` - Applied constant temperature to node ${r+1} (element ${n+1}, local node ${i+1})`),e[r]=1,t[r]=s}))}}))}})):"2D"===this.meshDimension&&Object.keys(this.boundaryConditions).forEach((n=>{if("constantTemp"===this.boundaryConditions[n][0]){const s=this.boundaryConditions[n][1];o(`Boundary ${n}: Applying constant temperature of ${s} K (Dirichlet condition)`),this.boundaryElements[n].forEach((([n,i])=>{if("linear"===this.elementOrder){({0:[0,2],1:[0,1],2:[1,3],3:[2,3]})[i].forEach((i=>{const r=this.nop[n][i]-1;o(` - Applied constant temperature to node ${r+1} (element ${n+1}, local node ${i+1})`),e[r]=1,t[r]=s}))}else if("quadratic"===this.elementOrder){({0:[0,3,6],1:[0,1,2],2:[2,5,8],3:[6,7,8]})[i].forEach((i=>{const r=this.nop[n][i]-1;o(` - Applied constant temperature to node ${r+1} (element ${n+1}, local node ${i+1})`),e[r]=1,t[r]=s}))}}))}}))}imposeConvectionBoundaryConditions(e,t,n,s,i,r,a){let l=[],d=[];Object.keys(this.boundaryConditions).forEach((e=>{const t=this.boundaryConditions[e];"convection"===t[0]&&(l[e]=t[1],d[e]=t[2])})),"1D"===this.meshDimension?Object.keys(this.boundaryConditions).forEach((n=>{if("convection"===this.boundaryConditions[n][0]){const s=l[n],i=d[n];o(`Boundary ${n}: Applying convection with heat transfer coefficient h=${s} W/(m²·K) and external temperature T∞=${i} K`),this.boundaryElements[n].forEach((([n,r])=>{let a;"linear"===this.elementOrder?a=0===r?0:1:"quadratic"===this.elementOrder&&(a=0===r?0:2);const l=this.nop[n][a]-1;o(` - Applied convection boundary condition to node ${l+1} (element ${n+1}, local node ${a+1})`),e[l]+=-s*i,t[l][l]+=s}))}})):"2D"===this.meshDimension&&Object.keys(this.boundaryConditions).forEach((c=>{if("convection"===this.boundaryConditions[c][0]){const u=l[c],m=d[c];o(`Boundary ${c}: Applying convection with heat transfer coefficient h=${u} W/(m²·K) and external temperature T∞=${m} K`),this.boundaryElements[c].forEach((([l,d])=>{if("linear"===this.elementOrder){let c,h,f,p,g;0===d?(c=n[0],h=0,f=0,p=3,g=2):1===d?(c=0,h=n[0],f=0,p=2,g=1):2===d?(c=n[0],h=1,f=1,p=4,g=2):3===d&&(c=1,h=n[0],f=2,p=4,g=1);let b=a.getBasisFunctions(c,h),y=b.basisFunction,v=b.basisFunctionDerivKsi,E=b.basisFunctionDerivEta,C=0,D=0,M=0,F=0;const $=this.nop[l].length;for(let e=0;e<$;e++){const t=this.nop[l][e]-1;0===d||2===d?(C+=i[t]*v[e],D+=r[t]*v[e]):1!==d&&3!==d||(M+=i[t]*E[e],F+=r[t]*E[e])}let w;w=0===d||2===d?Math.sqrt(C**2+D**2):Math.sqrt(M**2+F**2);for(let n=f;n{const t=this.boundaryConditions[e];"convection"===t[0]&&(a[e]=t[1],l[e]=t[2])}));const d=this.nop[e].length,c=Array(d).fill().map((()=>Array(d).fill(0))),u=Array(d).fill(0);for(const m in this.boundaryElements)if("convection"===this.boundaryConditions[m]?.[0]){const h=a[m],f=l[m];o(`Boundary ${m}: Applying convection with heat transfer coefficient h=${h} W/(m²·K) and external temperature T∞=${f} K`);const p=this.boundaryElements[m].find((([t,n])=>t===e));if(p){const a=p[1];if("1D"===this.meshDimension){let t;"linear"===this.elementOrder?t=0===a?0:1:"quadratic"===this.elementOrder&&(t=0===a?0:2),o(` - Applied convection boundary condition to node ${t+1} (element ${e+1}, local node ${t+1})`),u[t]+=-h*f,c[t][t]+=h}else if("2D"===this.meshDimension)if("linear"===this.elementOrder){let o,l,m,p,g;0===a?(o=s[0],l=0,m=0,p=3,g=2):1===a?(o=0,l=s[0],m=0,p=2,g=1):2===a?(o=s[0],l=1,m=1,p=4,g=2):3===a&&(o=1,l=s[0],m=2,p=4,g=1);const b=r.getBasisFunctions(o,l),y=b.basisFunction,v=b.basisFunctionDerivKsi,E=b.basisFunctionDerivEta;let C,D=0,M=0,F=0,$=0;for(let o=0;oArray(a).fill(0))),m=Array(a).fill(0),h=Array(a),f=Array(a);for(let n=0;ne-1)),{detJacobian:f,basisFunctionDerivX:p,basisFunctionDerivY:g}=S({basisFunction:n,basisFunctionDerivKsi:s,basisFunctionDerivEta:c,nodesXCoordinates:l,nodesYCoordinates:d,localToGlobalMap:m,nodesPerElement:a});for(let n=0;n{if("constantValue"===this.boundaryConditions[n][0]){const s=this.boundaryConditions[n][1];o(`Boundary ${n}: Applying constant value of ${s} (Dirichlet condition)`),this.boundaryElements[n].forEach((([n,i])=>{if("linear"===this.elementOrder){({0:[0],1:[1]})[i].forEach((i=>{const r=this.nop[n][i]-1;o(` - Applied constant value to node ${r+1} (element ${n+1}, local node ${i+1})`),e[r]=s;for(let n=0;n{const r=this.nop[n][i]-1;o(` - Applied constant value to node ${r+1} (element ${n+1}, local node ${i+1})`),e[r]=s;for(let n=0;n{if("constantValue"===this.boundaryConditions[n][0]){const s=this.boundaryConditions[n][1];o(`Boundary ${n}: Applying constant value of ${s} (Dirichlet condition)`),this.boundaryElements[n].forEach((([n,i])=>{if("linear"===this.elementOrder){({0:[0,2],1:[0,1],2:[1,3],3:[2,3]})[i].forEach((i=>{const r=this.nop[n][i]-1;o(` - Applied constant value to node ${r+1} (element ${n+1}, local node ${i+1})`),e[r]=s;for(let n=0;n{const r=this.nop[n][i]-1;o(` - Applied constant value to node ${r+1} (element ${n+1}, local node ${i+1})`),e[r]=s;for(let n=0;n{if("constantValue"===this.boundaryConditions[n][0]){const s=this.boundaryConditions[n][1];o(`Boundary ${n}: Applying constant value of ${s} (Dirichlet condition)`),this.boundaryElements[n].forEach((([n,i])=>{if("linear"===this.elementOrder){({0:[0],1:[1]})[i].forEach((i=>{const r=this.nop[n][i]-1;o(` - Applied constant value to node ${r+1} (element ${n+1}, local node ${i+1})`),e[r]=1,t[r]=s}))}else if("quadratic"===this.elementOrder){({0:[0],2:[2]})[i].forEach((i=>{const r=this.nop[n][i]-1;o(` - Applied constant value to node ${r+1} (element ${n+1}, local node ${i+1})`),e[r]=1,t[r]=s}))}}))}})):"2D"===this.meshDimension&&Object.keys(this.boundaryConditions).forEach((n=>{if("constantValue"===this.boundaryConditions[n][0]){const s=this.boundaryConditions[n][1];o(`Boundary ${n}: Applying constant value of ${s} (Dirichlet condition)`),this.boundaryElements[n].forEach((([n,i])=>{if("linear"===this.elementOrder){({0:[0,2],1:[0,1],2:[1,3],3:[2,3]})[i].forEach((i=>{const r=this.nop[n][i]-1;o(` - Applied constant value to node ${r+1} (element ${n+1}, local node ${i+1})`),e[r]=1,t[r]=s}))}else if("quadratic"===this.elementOrder){({0:[0,3,6],1:[0,1,2],2:[2,5,8],3:[6,7,8]})[i].forEach((i=>{const r=this.nop[n][i]-1;o(` - Applied constant value to node ${r+1} (element ${n+1}, local node ${i+1})`),e[r]=1,t[r]=s}))}}))}}))}}function K(e,t,n,i){s("Starting front propagation matrix assembly...");let r=1-i+.01;o(`eikonalViscousTerm: ${r}`),o(`eikonalActivationFlag: ${i}`);const{nodesXCoordinates:a,nodesYCoordinates:l,nop:d,boundaryElements:c,totalElements:u,meshDimension:m,elementOrder:h}=e,f=X(e),{residualVector:p,jacobianMatrix:g,localToGlobalMap:b,basisFunctions:y,gaussPoints:v,gaussWeights:E,nodesPerElement:C}=f;for(let e=0;eArray(d).fill(0))),p=Array(d).fill(0),g=Array(d),b=Array(d);for(let n=0;nArray(e).fill(0))),J.nodeConstraintCode=Array(e).fill(0),J.boundaryValues=Array(e).fill(0),J.globalResidualVector=Array(e).fill(0),J.solutionVector=Array(e).fill(0),J.topologyData=Array(t).fill(0),J.lateralData=Array(t).fill(0),H.writeFlag=0,H.totalNodes=e,H.transformationFlag=0,H.nodesPerElement=Array(t).fill(0),H.determinant=1;const n=Math.max(e,2e3);H.globalSolutionVector=Array(n).fill(0),H.frontDataIndex=0,U.localJacobianMatrix=Array(e).fill().map((()=>Array(e).fill(0))),U.currentElementIndex=0;const o=function(e,t){const n=Math.max(Math.ceil(Math.sqrt(t))*e,2*e);return n*t}(e,t);z.frontValues=Array(o).fill(0),z.columnHeaders=Array(n).fill(0),z.pivotRow=Array(n).fill(0),z.pivotData=Array(o).fill(0)}(a.nodesPerElement,d),s("Solving system using frontal..."),console.time("systemSolving"),_=new x({meshDimension:t.meshDimension,elementOrder:t.elementOrder});for(let e=0;eArray(l).fill(0))),y=Array(a).fill(0),v=Array(a).fill(0),E=Array(a).fill(0),C=1;H.writeFlag++;let D=1,M=1;U.currentElementIndex=0;for(let e=0;el||F>l)return void i("Error: systemSize not large enough");for(let e=0;e0)for(let e=0;eM||U.currentElementIndexMath.abs(n)&&(n=r,t=s,e=i)}}}let s=Math.abs(m[e-1]);d=Math.abs(z.columnHeaders[t-1]);let a=s+d+y[s-1]+v[d-1];H.determinant=H.determinant*n*(-1)**a/Math.abs(n);for(let e=0;e=s&&y[e]--,e>=d&&v[e]--;if(Math.abs(n)<1e-10&&i(`Matrix singular or ill-conditioned, currentElementIndex=${U.currentElementIndex}, pivotGlobalRowIndex=${s}, pivotColumnGlobalIndex=${d}, pivotValue=${n}`),0===n)return;for(let t=0;t1)for(let n=0;n1&&0!==o)for(let e=0;e1)for(let e=0;e1||U.currentElementIndex=e.totalElements)return i(`Skipping out-of-range elementIndex=${s} (totalElements=${e.totalElements})`),!1;const{localJacobianMatrix:r,localResidualVector:a,ngl:l}=o({elementIndex:s,nop:J.nodalNumbering,meshData:e,basisFunctions:_,FEAData:t,solutionVector:H.currentSolutionVector,eikonalActivationFlag:H.eikonalActivationFlag});let d=Array(t.nodesPerElement).fill().map((()=>Array(t.nodesPerElement).fill(0))),c=Array(t.nodesPerElement).fill(0);if(o===W){let o=!1;for(const t in e.boundaryElements)if("convection"===n.boundaryConditions[t]?.[0]&&e.boundaryElements[t].some((([e,t])=>e===s))){o=!0;break}if(o){const{gaussPoints:o,gaussWeights:i}=t,r=n.imposeConvectionBoundaryConditionsFront(s,e.nodesXCoordinates,e.nodesYCoordinates,o,i,_);d=r.localJacobianMatrix,c=r.localResidualVector}}for(let e=0;e0)continue;let r=0;z.pivotRow[s-1]=0;for(let e=0;e100){i(`Solution not converged. Error norm: ${o}`);break}a++}return{solutionVector:d,converged:r,iterations:a,jacobianMatrix:c,residualVector:u}}class ne{constructor(e,t,n,o,s,i,r,a){this.boundaryConditions=e,this.boundaryElements=t,this.nop=n,this.meshDimension=o,this.elementOrder=s,this.totalNodesVelocity=i,this.totalNodesPressure=r,this.q2ToPressureMap=a}imposeDirichletBoundaryConditions(e,t){const n=e.length;let s=!1;if("2D"===this.meshDimension&&(Object.keys(this.boundaryConditions).forEach((i=>{const r=this.boundaryConditions[i][0];if("stressFree"===r)s=!0,o(`Boundary ${i}: Applying stress-free condition (natural BC)`);else if("constantVelocity"===r){const s=this.boundaryConditions[i][1],r=this.boundaryConditions[i][2];o(`Boundary ${i}: Applying constant velocity condition (u=${s}, v=${r})`),this.boundaryElements[i].forEach((([i,a])=>{if("quadratic"===this.elementOrder){({0:[0,3,6],1:[0,1,2],2:[2,5,8],3:[6,7,8]})[a].forEach((a=>{const l=this.nop[i][a]-1,d=l,c=this.totalNodesVelocity+l;o(` - Applied velocity Dirichlet to node ${l+1} (element ${i+1}, local node ${a+1})`),e[d]=s;for(let e=0;e{const l=this.nop[i][a]-1,d=l,c=this.totalNodesVelocity+l;o(` - Applied velocity Dirichlet to node ${l+1} (element ${i+1}, local node ${a+1})`),e[d]=s;for(let e=0;e0&&(i.initialSolution=[...r]);const e=te(K,i);t=e.jacobianMatrix,n=e.residualVector,r=e.solutionVector,o+=1/s}}else if("generalFormPDEScript"===this.solverConfig)if("frontal"===this.solverMethod)i("Frontal solver is not yet supported for generalFormPDEScript. Please use 'lusolve' or 'jacobi'.");else{({jacobianMatrix:t,residualVector:n}=function(e,t,n){s("Starting general form PDE matrix assembly...");const{nodesXCoordinates:o,nodesYCoordinates:r,nop:a,boundaryElements:l,totalElements:d,meshDimension:c,elementOrder:u}=e,{A:m,B:h,C:f,D:p}=n,g=X(e),{residualVector:b,jacobianMatrix:y,localToGlobalMap:v,basisFunctions:E,gaussPoints:C,gaussWeights:D,nodesPerElement:M}=g;if("1D"===c)for(let e=0;e{let t={nodesXCoordinates:[],nodesYCoordinates:[],nodalNumbering:{quadElements:[],triangleElements:[]},boundaryElements:[],boundaryConditions:[],boundaryNodePairs:{},gmshV:0,ascii:!1,fltBytes:"8",totalNodesX:0,totalNodesY:0,physicalPropMap:[],elementTypes:{}};const n={curves:{}};let s=(await e.text()).split("\n").map((e=>e.trim())).filter((e=>""!==e)),i="",r=0,a=0,l=0,d=0,c={numNodes:0},u=[],m=0,h=0,f=0,p=0,g={numElements:0},b=0,y={},v=null,E=null,C=0,D=0,M=0;for(;rparseInt(e,10)));n.curves[e]=s,D++,D===v.curves&&(E="surfaces")}else"surfaces"===E&&(M++,M===v.surfaces&&(E="volumes"));else if("nodes"===i){if(0===a)a=parseInt(o[0],10),l=parseInt(o[1],10),t.nodesXCoordinates=new Array(l).fill(0),t.nodesYCoordinates=new Array(l).fill(0);else if(dparseInt(e,10)));let s=g.tag;if(1===g.dim){const e=n.curves[g.tag];e&&e.length>0&&(s=e[0])}1===g.elementType||8===g.elementType?(y[s]||(y[s]=[]),y[s].push(e),t.boundaryNodePairs[s]||(t.boundaryNodePairs[s]=[]),t.boundaryNodePairs[s].push(e)):2===g.elementType?t.nodalNumbering.triangleElements.push(e):3!==g.elementType&&10!==g.elementType||t.nodalNumbering.quadElements.push(e),b++,b===g.numElements&&(p++,g={numElements:0})}r++}return t.physicalPropMap.forEach((e=>{if(1===e.dimension){const n=y[e.tag]||[];n.length>0&&t.boundaryConditions.push({name:e.name,tag:e.tag,nodes:n})}})),o(`Parsed boundary node pairs by physical tag: ${JSON.stringify(t.boundaryNodePairs)}`),t};function ie(e,t,n,o){console.time("plottingTime");const{nodesXCoordinates:s,nodesYCoordinates:i}=t.nodesCoordinates,r=t.solutionVector,a=e.solverConfig,l=e.meshConfig.meshDimension;if(V(e.meshConfig),"1D"===l&&"line"===n){let e;e=r.length>0&&Array.isArray(r[0])?r.map((e=>e[0])):r,Array.from(s);let t={x:s,y:e,mode:"lines",type:"scatter",line:{color:"rgb(219, 64, 82)",width:2},name:"Solution"},n=Math.min(window.innerWidth,700),i={title:`line plot - ${a}`,width:Math.min(n,600),height:300,xaxis:{title:"x"},yaxis:{title:"Solution"},margin:{l:50,r:50,t:50,b:50}};Plotly.newPlot(o,[t],i,{responsive:!0}),console.timeEnd("plottingTime")}else if("2D"===l&&"contour"===n){let e;e=Array.isArray(r[0])?r.map((e=>e[0])):r;let t=Math.min(window.innerWidth,700),l=Math.min(...s),d=Math.max(...s),c=Math.min(...i),u=(Math.max(...i)-c)/(d-l),m=Math.min(t,600),h={title:`${n} plot - ${a}`,width:m,height:m*u,xaxis:{title:"x"},yaxis:{title:"y",scaleanchor:"x",scaleratio:1},margin:{l:50,r:50,t:50,b:50},hovermode:"closest"},f={x:s,y:i,z:e,type:"contour",line:{smoothing:.85},contours:{coloring:"heatmap",showlabels:!1},colorbar:{title:"Solution"},name:"Solution Field"};Plotly.newPlot(o,[f],h,{responsive:!0}),console.timeEnd("plottingTime")}}function re(e,t,n,o){console.time("plottingTime");const{nodesXCoordinates:s,nodesYCoordinates:i}=t.nodesCoordinates,r=e.meshConfig.meshDimension,a=V(e.meshConfig),l=new x({meshDimension:e.meshConfig.meshDimension,elementOrder:e.meshConfig.elementOrder});if("1D"===r&&"line"===n);else if("2D"===r&&"contour"===n){const r=[],d=[],c=Math.max(...s)-Math.min(...s),u=Math.max(...i)-Math.min(...i),m=50,h=Math.round(c*m),f=Math.round(u*m),p=c/(h-1),g=u/(f-1);let b=[];r[0]=Math.min(...s),d[0]=Math.min(...i);for(let e=1;ee[0])):a;let u=0;for(let e=0;et!=l>t&&e<(a-i)*(t-r)/(l-r)+i&&(o=!o)}return o}let ce=null;async function ue(){if(ce)return ce;await import("@kitware/vtk.js/Rendering/Profiles/Geometry.js");const[{default:e},{default:t},{default:n},{default:o},{default:s},{default:i},{default:r},{default:a},{default:l},{default:d}]=await Promise.all([import("@kitware/vtk.js/Rendering/Core/Actor.js"),import("@kitware/vtk.js/Rendering/Core/ColorTransferFunction.js"),import("@kitware/vtk.js/Rendering/Core/ColorTransferFunction/ColorMaps.js"),import("@kitware/vtk.js/Common/Core/DataArray.js"),import("@kitware/vtk.js/Common/DataModel/ImageData.js"),import("@kitware/vtk.js/Filters/General/ImageMarchingSquares.js"),import("@kitware/vtk.js/Rendering/Misc/GenericRenderWindow.js"),import("@kitware/vtk.js/Rendering/Core/Mapper.js"),import("@kitware/vtk.js/Common/DataModel/PolyData.js"),import("@kitware/vtk.js/Rendering/Core/ScalarBarActor.js")]);return ce={vtkActor:e,vtkColorTransferFunction:t,vtkColorMaps:n,vtkDataArray:o,vtkImageData:s,vtkImageMarchingSquares:i,vtkGenericRenderWindow:r,vtkMapper:a,vtkPolyData:l,vtkScalarBarActor:d},ce}const me=new Map;function he(e={}){return{presetName:e.presetName??"Cool to Warm",reverse:e.reverse??!1,showScalarBar:e.showScalarBar??!0,scalarBarTitle:e.scalarBarTitle??"Solution"}}function fe(e={}){return{enabled:e.enabled??!0,numberOfContours:e.numberOfContours??12,color:e.color??[.15,.15,.15],lineWidth:e.lineWidth??1}}async function pe(e,t,n,o,s={}){console.time("plottingTime");const i=e.meshConfig.meshDimension,r=V(e.meshConfig),a=await be(e,t,r,{mode:"1D"===i&&"line"===n?"line":"surface"});await Ee(a,o,e.solverConfig,n,s),console.timeEnd("plottingTime")}async function ge(e,t,n,o,s={}){if(console.time("plottingTime"),"2D"!==e.meshConfig.meshDimension||"contour"!==n){const i=V(e.meshConfig),r=await be(e,t,i,{mode:"1D"===e.meshConfig.meshDimension&&"line"===n?"line":"surface"});return await Ee(r,o,e.solverConfig,n,s),void console.timeEnd("plottingTime")}const i=V(e.meshConfig),r=await async function(e,t,n){const o=await async function(e,t,n){const{nodesXCoordinates:o,nodesYCoordinates:s}=t.nodesCoordinates,i=new x({meshDimension:e.meshConfig.meshDimension,elementOrder:e.meshConfig.elementOrder});let r=o[0],a=o[0],l=s[0],d=s[0];for(let e=1;ea&&(a=t),nd&&(d=n)}const c=a-r,u=d-l,m=50,h=Math.max(2,Math.round(c*m)),f=Math.max(2,Math.round(u*m)),p=c/(h-1),g=u/(f-1),b=h*f,y=new Float32Array(b),v=new Float32Array(b),E=new Float32Array(b);E.fill(Number.NaN);const C=new Uint8Array(b),D=B(n),{nodeNeighbors:M,neighborCount:F}=Y(n);let $=0;for(let o=0;o',''," ",` `,' ',` ${Array.from(e.scalars).join(" ")}`," "," ",` ${Array.from(e.points).join(" ")}`," ",` <${i}>`,` ${t.join(" ")}`,` ${n.join(" ")}`,` ${i}>`," "," ",""].join("\n")}(await be(e,t,n,o))}function ve(e){const{nodesXCoordinates:t,nodesYCoordinates:n}=e.nodesCoordinates,o=Me(e.solutionVector,t.length),s=new Float32Array(3*t.length);for(let e=0;en&&(n=s))}return Number.isFinite(t)&&Number.isFinite(n)?t===n?[t-1,n+1]:[t,n]:[0,1]}(e.scalars);g.setScalarRange(b[0],b[1]);const y=s.colorScale??he({}),v=a.newInstance(),E=function(e,t){if(!t||!e?.RGBPoints?.length)return e;const n=e.RGBPoints,o=n[0],s=n[n.length-4],i=[];for(let e=n.length-4;e>=0;e-=4)i.push(s-(n[e]-o),n[e+1],n[e+2],n[e+3]);return{...e,RGBPoints:i}}(l.getPresetByName(y.presetName)??l.getPresetByName("Cool to Warm"),y.reverse);v.applyColorMap(E),v.setMappingRange(b[0],b[1]),v.updateRange(),g.setLookupTable(v);const C=r.newInstance();if(C.setMapper(g),"line"===e.mode&&C.getProperty().setLineWidth(3),f.addActor(C),y.showScalarBar){const e=u.newInstance({drawNanAnnotation:!1,generateTicks:e=>{const t=e.getLastTickBounds();if(!t||t.length<2)return;const[n,o]=t,s=(o-n)/4,i=Array.from({length:5},((e,t)=>n+t*s));e.setTicks(i),e.setTickStrings(i.map((e=>{const t=Math.abs(e);return 0===t?"0":t>=.01&&t<1e4?parseFloat(e.toPrecision(4)).toString():e.toExponential(2)})))}});e.setTickTextStyle({fontColor:"black"}),e.setAxisTextStyle({fontColor:"black"}),e.setAxisLabel(y.scalarBarTitle),e.setScalarsToColors(v),f.addActor2D(e)}const D=fe(s.contourLines??{enabled:!1});D.enabled&&"line"!==e.mode&&await async function(e,t,n,o){const s=t.metadata?.interpolationGrid;if(!s)return;const{vtkActor:i,vtkDataArray:r,vtkImageData:a,vtkImageMarchingSquares:l,vtkMapper:d}=await ue(),c=a.newInstance();c.setDimensions(s.nx,s.ny,1),c.setOrigin(s.origin[0],s.origin[1],0),c.setSpacing(s.spacing[0],s.spacing[1],1);const u=r.newInstance({name:"solution",numberOfComponents:1,values:s.imageScalars});c.getPointData().setScalars(u);const m=l.newInstance({slicingMode:2,slice:0,mergePoints:!0});m.setInputData(c);const h=Math.max(2,o.numberOfContours),f=n[0],p=n[1],g=(p-f)/(h-1),b=[];for(let e=0;e0&&("line"===o?r.getLines().setData(n):r.getPolys().setData(n));const a=i.newInstance({name:"solution",numberOfComponents:1,values:t});return r.getPointData().setScalars(a),r}function De(e,t){const n=new Float32Array(3*e.length);for(let o=0;oe-1)),n=t.length;return 2===n||3===n?t:4===n?[t[0],t[2],t[3],t[1]]:6===n?[t[0],t[2],t[5]]:8===n?[t[0],t[4],t[6],t[2]]:9===n?[t[0],t[6],t[8],t[2]]:t.slice(0,Math.min(4,t.length))}function $e(e){let t=0,n=0;for(;ne[0])):a;let u=0;for(let e=0;et!=l>t&&e<(a-i)*(t-r)/(l-r)+i&&(o=!o)}return o}class ke{constructor(){this.worker=null,this.feaWorker=null,this.isReady=!1,this._initWorker()}async _initWorker(){try{this.worker=new Worker(new URL("./wrapperScript.js",import.meta.url),{type:"module"}),this.worker.onerror=e=>{console.error("FEAScriptWorker: Worker error:",e)};const e=p(this.worker);this.feaWorker=await new e,this.isReady=!0}catch(e){throw console.error("Failed to initialize worker",e),e}}async _ensureReady(){return this.isReady?Promise.resolve():new Promise(((e,t)=>{let n=0;const o=()=>{n++,this.isReady?e():n>=50?t(new Error("Timeout waiting for worker to be ready")):setTimeout(o,1e3)};o()}))}async setModelConfig(e){return await this._ensureReady(),this.feaWorker.setModelConfig(e)}async setMeshConfig(e){return await this._ensureReady(),this.feaWorker.setMeshConfig(e)}async addBoundaryCondition(e,t){return await this._ensureReady(),this.feaWorker.addBoundaryCondition(e,t)}async setSolverMethod(e){return await this._ensureReady(),this.feaWorker.setSolverMethod(e)}async solve(){await this._ensureReady(),performance.now();const e=await this.feaWorker.solve();return performance.now(),e}async getModelInfo(){return await this._ensureReady(),this.feaWorker.getModelInfo()}async ping(){return await this._ensureReady(),this.feaWorker.ping()}terminate(){this.worker&&(this.worker.terminate(),this.worker=null,this.feaWorker=null,this.isReady=!1)}}const Pe="0.3.0 (RC)";export{oe as FEAScriptModel,ke as FEAScriptWorker,he as createColorScale,fe as createContourLineOptions,se as importGmshQuadTri,n as logSystem,re as plotInterpolatedSolution,ge as plotInterpolatedSolutionVtk,ie as plotSolution,pe as plotSolutionVtk,Pe as printVersion,ve as transformSolverOutputToMLBuffers,ye as transformSolverOutputToVTP,be as transformSolverOutputToVtkData};
//# sourceMappingURL=feascript.esm.js.map
diff --git a/dist/feascript.esm.js.map b/dist/feascript.esm.js.map
index d14ea29..1eb3e89 100644
--- a/dist/feascript.esm.js.map
+++ b/dist/feascript.esm.js.map
@@ -1 +1 @@
-{"version":3,"file":"feascript.esm.js","sources":["../src/methods/euclideanNormScript.js","../src/utilities/loggingScript.js","../src/vendor/comlink.mjs","../src/methods/linearSystemSolverScript.js","../src/methods/jacobiSolverScript.js","../src/mesh/basisFunctionsScript.js","../src/mesh/meshGenerationScript.js","../src/methods/numericalIntegrationScript.js","../src/mesh/meshUtilsScript.js","../src/models/thermalBoundaryConditionsScript.js","../src/models/heatConductionScript.js","../src/models/genericBoundaryConditionsScript.js","../src/models/frontPropagationScript.js","../src/methods/frontalSolverScript.js","../src/methods/newtonRaphsonScript.js","../src/models/flowBoundaryConditions.js","../src/FEAScript.js","../src/models/generalFormPDEScript.js","../src/models/creepingFlowScript.js","../src/readers/gmshReaderScript.js","../src/visualization/plotSolutionScript.js","../src/workers/workerScript.js","../src/index.js"],"sourcesContent":["/**\n * ════════════════════════════════════════════════════════════════\n * FEAScript Core Library\n * Lightweight Finite Element Simulation in JavaScript\n * Version: 0.2.0 | https://feascript.com\n * MIT License © 2023–2026 FEAScript\n * ════════════════════════════════════════════════════════════════\n */\n\n/**\n * Function to calculate the Euclidean norm of a vector\n * @param {array} vector - The input vector\n * @returns {number} The Euclidean norm of the vector\n */\nexport function euclideanNorm(vector) {\n let norm = 0;\n for (let i = 0; i < vector.length; i++) {\n norm += vector[i] * vector[i];\n }\n norm = Math.sqrt(norm);\n return norm;\n}\n","/**\n * ════════════════════════════════════════════════════════════════\n * FEAScript Core Library\n * Lightweight Finite Element Simulation in JavaScript\n * Version: 0.2.0 | https://feascript.com\n * MIT License © 2023–2026 FEAScript\n * ════════════════════════════════════════════════════════════════\n */\n\n// Global logging level\nlet currentLogLevel = \"basic\";\n\n/**\n * Function to set the logging system level\n * @param {string} level - Logging level (basic, debug)\n */\nexport function logSystem(level) {\n if (level !== \"basic\" && level !== \"debug\") {\n console.log(\n \"%c[WARN] Invalid log level: \" + level + \". Using basic instead.\",\n \"color: #FFC107; font-weight: bold;\"\n ); // Yellow for warnings\n currentLogLevel = \"basic\";\n } else {\n currentLogLevel = level;\n basicLog(`Log level set to: ${level}`);\n }\n}\n\n/**\n * Function to log debug messages - only logs if level is 'debug'\n * @param {string} message - Message to log\n */\nexport function debugLog(message) {\n if (currentLogLevel === \"debug\") {\n console.log(\"%c[DEBUG] \" + message, \"color: #2196F3; font-weight: bold;\");\n }\n}\n\n/**\n * Function to log basic information - always logs\n * @param {string} message - Message to log\n */\nexport function basicLog(message) {\n console.log(\"%c[INFO] \" + message, \"color: #4CAF50; font-weight: bold;\");\n}\n\n/**\n * Function to log error messages\n * @param {string} message - Message to log\n */\nexport function errorLog(message) {\n console.log(\"%c[ERROR] \" + message, \"color: #F44336; font-weight: bold;\");\n}\n\n/**\n * Function to log warning messages\n * @param {string} message - Message to log\n */\nexport function warnLog(message) {\n console.log(\"%c[WARN] \" + message, \"color: #FF9800; font-weight: bold;\");\n}\n\n/**\n * Function to handle version information and fetch the latest update date and release from GitHub\n */\nexport async function printVersionInformation() {\n basicLog(\"Fetching latest FEAScript version information...\");\n try {\n const commitResponse = await fetch(\"https://api.github.com/repos/FEAScript/FEAScript/commits/main\");\n const commitData = await commitResponse.json();\n const latestCommitDate = new Date(commitData.commit.committer.date).toLocaleString();\n basicLog(`Latest FEAScript update: ${latestCommitDate}`);\n return latestCommitDate;\n } catch (error) {\n errorLog(\"Failed to fetch version information: \" + error);\n return \"Version information unavailable\";\n }\n}\n","/**\n * @license\n * Copyright 2019 Google LLC\n * SPDX-License-Identifier: Apache-2.0\n */\nconst proxyMarker = Symbol(\"Comlink.proxy\");\nconst createEndpoint = Symbol(\"Comlink.endpoint\");\nconst releaseProxy = Symbol(\"Comlink.releaseProxy\");\nconst finalizer = Symbol(\"Comlink.finalizer\");\nconst throwMarker = Symbol(\"Comlink.thrown\");\nconst isObject = (val) => (typeof val === \"object\" && val !== null) || typeof val === \"function\";\n/**\n * Internal transfer handle to handle objects marked to proxy.\n */\nconst proxyTransferHandler = {\n canHandle: (val) => isObject(val) && val[proxyMarker],\n serialize(obj) {\n const { port1, port2 } = new MessageChannel();\n expose(obj, port1);\n return [port2, [port2]];\n },\n deserialize(port) {\n port.start();\n return wrap(port);\n },\n};\n/**\n * Internal transfer handler to handle thrown exceptions.\n */\nconst throwTransferHandler = {\n canHandle: (value) => isObject(value) && throwMarker in value,\n serialize({ value }) {\n let serialized;\n if (value instanceof Error) {\n serialized = {\n isError: true,\n value: {\n message: value.message,\n name: value.name,\n stack: value.stack,\n },\n };\n }\n else {\n serialized = { isError: false, value };\n }\n return [serialized, []];\n },\n deserialize(serialized) {\n if (serialized.isError) {\n throw Object.assign(new Error(serialized.value.message), serialized.value);\n }\n throw serialized.value;\n },\n};\n/**\n * Allows customizing the serialization of certain values.\n */\nconst transferHandlers = new Map([\n [\"proxy\", proxyTransferHandler],\n [\"throw\", throwTransferHandler],\n]);\nfunction isAllowedOrigin(allowedOrigins, origin) {\n for (const allowedOrigin of allowedOrigins) {\n if (origin === allowedOrigin || allowedOrigin === \"*\") {\n return true;\n }\n if (allowedOrigin instanceof RegExp && allowedOrigin.test(origin)) {\n return true;\n }\n }\n return false;\n}\nfunction expose(obj, ep = globalThis, allowedOrigins = [\"*\"]) {\n ep.addEventListener(\"message\", function callback(ev) {\n if (!ev || !ev.data) {\n return;\n }\n if (!isAllowedOrigin(allowedOrigins, ev.origin)) {\n console.warn(`Invalid origin '${ev.origin}' for comlink proxy`);\n return;\n }\n const { id, type, path } = Object.assign({ path: [] }, ev.data);\n const argumentList = (ev.data.argumentList || []).map(fromWireValue);\n let returnValue;\n try {\n const parent = path.slice(0, -1).reduce((obj, prop) => obj[prop], obj);\n const rawValue = path.reduce((obj, prop) => obj[prop], obj);\n switch (type) {\n case \"GET\" /* MessageType.GET */:\n {\n returnValue = rawValue;\n }\n break;\n case \"SET\" /* MessageType.SET */:\n {\n parent[path.slice(-1)[0]] = fromWireValue(ev.data.value);\n returnValue = true;\n }\n break;\n case \"APPLY\" /* MessageType.APPLY */:\n {\n returnValue = rawValue.apply(parent, argumentList);\n }\n break;\n case \"CONSTRUCT\" /* MessageType.CONSTRUCT */:\n {\n const value = new rawValue(...argumentList);\n returnValue = proxy(value);\n }\n break;\n case \"ENDPOINT\" /* MessageType.ENDPOINT */:\n {\n const { port1, port2 } = new MessageChannel();\n expose(obj, port2);\n returnValue = transfer(port1, [port1]);\n }\n break;\n case \"RELEASE\" /* MessageType.RELEASE */:\n {\n returnValue = undefined;\n }\n break;\n default:\n return;\n }\n }\n catch (value) {\n returnValue = { value, [throwMarker]: 0 };\n }\n Promise.resolve(returnValue)\n .catch((value) => {\n return { value, [throwMarker]: 0 };\n })\n .then((returnValue) => {\n const [wireValue, transferables] = toWireValue(returnValue);\n ep.postMessage(Object.assign(Object.assign({}, wireValue), { id }), transferables);\n if (type === \"RELEASE\" /* MessageType.RELEASE */) {\n // detach and deactive after sending release response above.\n ep.removeEventListener(\"message\", callback);\n closeEndPoint(ep);\n if (finalizer in obj && typeof obj[finalizer] === \"function\") {\n obj[finalizer]();\n }\n }\n })\n .catch((error) => {\n // Send Serialization Error To Caller\n const [wireValue, transferables] = toWireValue({\n value: new TypeError(\"Unserializable return value\"),\n [throwMarker]: 0,\n });\n ep.postMessage(Object.assign(Object.assign({}, wireValue), { id }), transferables);\n });\n });\n if (ep.start) {\n ep.start();\n }\n}\nfunction isMessagePort(endpoint) {\n return endpoint.constructor.name === \"MessagePort\";\n}\nfunction closeEndPoint(endpoint) {\n if (isMessagePort(endpoint))\n endpoint.close();\n}\nfunction wrap(ep, target) {\n const pendingListeners = new Map();\n ep.addEventListener(\"message\", function handleMessage(ev) {\n const { data } = ev;\n if (!data || !data.id) {\n return;\n }\n const resolver = pendingListeners.get(data.id);\n if (!resolver) {\n return;\n }\n try {\n resolver(data);\n }\n finally {\n pendingListeners.delete(data.id);\n }\n });\n return createProxy(ep, pendingListeners, [], target);\n}\nfunction throwIfProxyReleased(isReleased) {\n if (isReleased) {\n throw new Error(\"Proxy has been released and is not useable\");\n }\n}\nfunction releaseEndpoint(ep) {\n return requestResponseMessage(ep, new Map(), {\n type: \"RELEASE\" /* MessageType.RELEASE */,\n }).then(() => {\n closeEndPoint(ep);\n });\n}\nconst proxyCounter = new WeakMap();\nconst proxyFinalizers = \"FinalizationRegistry\" in globalThis &&\n new FinalizationRegistry((ep) => {\n const newCount = (proxyCounter.get(ep) || 0) - 1;\n proxyCounter.set(ep, newCount);\n if (newCount === 0) {\n releaseEndpoint(ep);\n }\n });\nfunction registerProxy(proxy, ep) {\n const newCount = (proxyCounter.get(ep) || 0) + 1;\n proxyCounter.set(ep, newCount);\n if (proxyFinalizers) {\n proxyFinalizers.register(proxy, ep, proxy);\n }\n}\nfunction unregisterProxy(proxy) {\n if (proxyFinalizers) {\n proxyFinalizers.unregister(proxy);\n }\n}\nfunction createProxy(ep, pendingListeners, path = [], target = function () { }) {\n let isProxyReleased = false;\n const proxy = new Proxy(target, {\n get(_target, prop) {\n throwIfProxyReleased(isProxyReleased);\n if (prop === releaseProxy) {\n return () => {\n unregisterProxy(proxy);\n releaseEndpoint(ep);\n pendingListeners.clear();\n isProxyReleased = true;\n };\n }\n if (prop === \"then\") {\n if (path.length === 0) {\n return { then: () => proxy };\n }\n const r = requestResponseMessage(ep, pendingListeners, {\n type: \"GET\" /* MessageType.GET */,\n path: path.map((p) => p.toString()),\n }).then(fromWireValue);\n return r.then.bind(r);\n }\n return createProxy(ep, pendingListeners, [...path, prop]);\n },\n set(_target, prop, rawValue) {\n throwIfProxyReleased(isProxyReleased);\n // FIXME: ES6 Proxy Handler `set` methods are supposed to return a\n // boolean. To show good will, we return true asynchronously ¯\\_(ツ)_/¯\n const [value, transferables] = toWireValue(rawValue);\n return requestResponseMessage(ep, pendingListeners, {\n type: \"SET\" /* MessageType.SET */,\n path: [...path, prop].map((p) => p.toString()),\n value,\n }, transferables).then(fromWireValue);\n },\n apply(_target, _thisArg, rawArgumentList) {\n throwIfProxyReleased(isProxyReleased);\n const last = path[path.length - 1];\n if (last === createEndpoint) {\n return requestResponseMessage(ep, pendingListeners, {\n type: \"ENDPOINT\" /* MessageType.ENDPOINT */,\n }).then(fromWireValue);\n }\n // We just pretend that `bind()` didn’t happen.\n if (last === \"bind\") {\n return createProxy(ep, pendingListeners, path.slice(0, -1));\n }\n const [argumentList, transferables] = processArguments(rawArgumentList);\n return requestResponseMessage(ep, pendingListeners, {\n type: \"APPLY\" /* MessageType.APPLY */,\n path: path.map((p) => p.toString()),\n argumentList,\n }, transferables).then(fromWireValue);\n },\n construct(_target, rawArgumentList) {\n throwIfProxyReleased(isProxyReleased);\n const [argumentList, transferables] = processArguments(rawArgumentList);\n return requestResponseMessage(ep, pendingListeners, {\n type: \"CONSTRUCT\" /* MessageType.CONSTRUCT */,\n path: path.map((p) => p.toString()),\n argumentList,\n }, transferables).then(fromWireValue);\n },\n });\n registerProxy(proxy, ep);\n return proxy;\n}\nfunction myFlat(arr) {\n return Array.prototype.concat.apply([], arr);\n}\nfunction processArguments(argumentList) {\n const processed = argumentList.map(toWireValue);\n return [processed.map((v) => v[0]), myFlat(processed.map((v) => v[1]))];\n}\nconst transferCache = new WeakMap();\nfunction transfer(obj, transfers) {\n transferCache.set(obj, transfers);\n return obj;\n}\nfunction proxy(obj) {\n return Object.assign(obj, { [proxyMarker]: true });\n}\nfunction windowEndpoint(w, context = globalThis, targetOrigin = \"*\") {\n return {\n postMessage: (msg, transferables) => w.postMessage(msg, targetOrigin, transferables),\n addEventListener: context.addEventListener.bind(context),\n removeEventListener: context.removeEventListener.bind(context),\n };\n}\nfunction toWireValue(value) {\n for (const [name, handler] of transferHandlers) {\n if (handler.canHandle(value)) {\n const [serializedValue, transferables] = handler.serialize(value);\n return [\n {\n type: \"HANDLER\" /* WireValueType.HANDLER */,\n name,\n value: serializedValue,\n },\n transferables,\n ];\n }\n }\n return [\n {\n type: \"RAW\" /* WireValueType.RAW */,\n value,\n },\n transferCache.get(value) || [],\n ];\n}\nfunction fromWireValue(value) {\n switch (value.type) {\n case \"HANDLER\" /* WireValueType.HANDLER */:\n return transferHandlers.get(value.name).deserialize(value.value);\n case \"RAW\" /* WireValueType.RAW */:\n return value.value;\n }\n}\nfunction requestResponseMessage(ep, pendingListeners, msg, transfers) {\n return new Promise((resolve) => {\n const id = generateUUID();\n pendingListeners.set(id, resolve);\n if (ep.start) {\n ep.start();\n }\n ep.postMessage(Object.assign({ id }, msg), transfers);\n });\n}\nfunction generateUUID() {\n return new Array(4)\n .fill(0)\n .map(() => Math.floor(Math.random() * Number.MAX_SAFE_INTEGER).toString(16))\n .join(\"-\");\n}\n\nexport { createEndpoint, expose, finalizer, proxy, proxyMarker, releaseProxy, transfer, transferHandlers, windowEndpoint, wrap };\n//# sourceMappingURL=comlink.mjs.map\n","/**\n * ════════════════════════════════════════════════════════════════\n * FEAScript Core Library\n * Lightweight Finite Element Simulation in JavaScript\n * Version: 0.2.0 | https://feascript.com\n * MIT License © 2023–2026 FEAScript\n * ════════════════════════════════════════════════════════════════\n */\n\n// Internal imports\nimport { jacobiSolver } from \"./jacobiSolverScript.js\";\nimport { basicLog, debugLog, errorLog } from \"../utilities/loggingScript.js\";\nimport * as Comlink from \"../vendor/comlink.mjs\";\n\n/**\n * Function to solve a system of linear equations using different solver methods\n * @param {string} solverMethod - The solver method to use (\"lusolve\" or \"jacobi\")\n * @param {Array} jacobianMatrix - The coefficient matrix\n * @param {Array} residualVector - The right-hand side vector\n * @param {object} [options] - Optional parameters for the solver, such as `maxIterations` and `tolerance`\n * @returns {object} An object containing:\n * - solutionVector: The solution vector\n * - converged: Boolean indicating whether the method converged (for iterative methods)\n * - iterations: Number of iterations performed (for iterative methods)\n */\nexport function solveLinearSystem(solverMethod, jacobianMatrix, residualVector, options = {}) {\n\n // Extract options\n const { maxIterations = 10000, tolerance = 1e-4 } = options;\n\n let solutionVector = [];\n let converged = true;\n let iterations = 0;\n\n // Solve the linear system based on the specified solver method\n basicLog(`Solving system using ${solverMethod}...`);\n console.time(\"systemSolving\");\n\n if (solverMethod === \"lusolve\") {\n // Use LU decomposition method\n const jacobianMatrixSparse = math.sparse(jacobianMatrix);\n const luFactorization = math.slu(jacobianMatrixSparse, 1, 1); // order=1, threshold=1 for pivoting\n let solutionMatrix = math.lusolve(luFactorization, residualVector);\n solutionVector = math.squeeze(solutionMatrix).valueOf();\n //solutionVector = math.lusolve(jacobianMatrix, residualVector); // In the case of a dense matrix\n } else if (solverMethod === \"jacobi\") {\n // Use Jacobi method\n const initialGuess = new Array(residualVector.length).fill(0);\n const jacobiSolverResult = jacobiSolver(jacobianMatrix, residualVector, initialGuess, {\n maxIterations,\n tolerance,\n });\n\n // Log convergence information\n if (jacobiSolverResult.converged) {\n debugLog(`Jacobi method converged in ${jacobiSolverResult.iterations} iterations`);\n } else {\n errorLog(`Jacobi method did not converge after ${jacobiSolverResult.iterations} iterations`);\n }\n\n solutionVector = jacobiSolverResult.solutionVector;\n converged = jacobiSolverResult.converged;\n iterations = jacobiSolverResult.iterations;\n } else {\n errorLog(`Unknown solver method: ${solverMethod}`);\n }\n\n console.timeEnd(\"systemSolving\");\n basicLog(\"System solved successfully\");\n\n return { solutionVector, converged, iterations };\n}\n\n// Helper to lazily create a default WebGPU compute engine (Comlink + worker)\nasync function createDefaultComputeEngine() {\n const worker = new Worker(new URL(\"../workers/webgpuWorkerScript.js\", import.meta.url), {\n type: \"module\",\n });\n const computeEngine = Comlink.wrap(worker);\n await computeEngine.initialize();\n return { computeEngine, worker };\n}\n\n/**\n * Function to solve asynchronously a system of linear equations using different solver methods\n * @param {string} solverMethod - The solver method to use (e.g., \"jacobi-gpu\")\n * @param {array} jacobianMatrix - The coefficient matrix\n * @param {array} residualVector - The right-hand side vector\n * @param {object} [options] - Optional parameters for the solver, such as `maxIterations` and `tolerance`\n * @returns {Promise