@@ -40,7 +40,23 @@ export class ThermalBoundaryConditions {
4040 const tempValue = this . boundaryConditions [ boundaryKey ] [ 1 ] ;
4141 this . boundaryElements [ boundaryKey ] . forEach ( ( [ elementIndex , side ] ) => {
4242 if ( this . elementOrder === "linear" ) {
43- console . log ( "Unsupported dimension or element order" ) ;
43+ const boundarySides = {
44+ 0 : [ 0 , 2 ] , // Nodes at the bottom side of the reference element
45+ 1 : [ 0 , 1 ] , // Nodes at the left side of the reference element
46+ 2 : [ 1 , 3 ] , // Nodes at the top side of the reference element
47+ 3 : [ 2 , 3 ] , // Nodes at the right side of the reference element
48+ } ;
49+ boundarySides [ boundaryKey ] . forEach ( ( nodeIndex ) => {
50+ const globalNodeIndex = this . nop [ elementIndex ] [ nodeIndex ] - 1 ;
51+ // Set the residual vector to the ConstantTemp value
52+ residualVector [ globalNodeIndex ] = tempValue ;
53+ // Set the Jacobian matrix row to zero
54+ for ( let colIndex = 0 ; colIndex < residualVector . length ; colIndex ++ ) {
55+ jacobianMatrix [ globalNodeIndex ] [ colIndex ] = 0 ;
56+ }
57+ // Set the diagonal entry of the Jacobian matrix to one
58+ jacobianMatrix [ globalNodeIndex ] [ globalNodeIndex ] = 1 ;
59+ } ) ;
4460 } else if ( this . elementOrder === "quadratic" ) {
4561 const boundarySides = {
4662 0 : [ 0 , 3 , 6 ] , // Nodes at the bottom side of the reference element
@@ -96,7 +112,107 @@ export class ThermalBoundaryConditions {
96112 const extTemp = convectionExtTemp [ boundaryKey ] ;
97113 this . boundaryElements [ boundaryKey ] . forEach ( ( [ elementIndex , side ] ) => {
98114 if ( this . elementOrder === "linear" ) {
99- console . log ( "Unsupported dimension or element order" ) ;
115+ let gaussPoint1 , gaussPoint2 , firstNodeIndex , lastNodeIndex , nodeIncrement ;
116+ if ( side === 0 ) {
117+ // Nodes at the bottom side of the reference element
118+ gaussPoint1 = gaussPoints [ 0 ] ;
119+ gaussPoint2 = 0 ;
120+ firstNodeIndex = 0 ;
121+ lastNodeIndex = 2 ;
122+ nodeIncrement = 2 ;
123+ } else if ( side === 1 ) {
124+ // Nodes at the left side of the reference element
125+ gaussPoint1 = 0 ;
126+ gaussPoint2 = gaussPoints [ 0 ] ;
127+ firstNodeIndex = 0 ;
128+ lastNodeIndex = 1 ;
129+ nodeIncrement = 1 ;
130+ } else if ( side === 2 ) {
131+ // Nodes at the top side of the reference element
132+ gaussPoint1 = gaussPoints [ 0 ] ;
133+ gaussPoint2 = 1 ;
134+ firstNodeIndex = 1 ;
135+ lastNodeIndex = 3 ;
136+ nodeIncrement = 2 ;
137+ } else if ( side === 3 ) {
138+ // Nodes at the right side of the reference element
139+ gaussPoint1 = 1 ;
140+ gaussPoint2 = gaussPoints [ 0 ] ;
141+ firstNodeIndex = 2 ;
142+ lastNodeIndex = 3 ;
143+ nodeIncrement = 1 ;
144+ }
145+ let basisFunctionsAndDerivatives = basisFunctionsData . getBasisFunctions (
146+ gaussPoint1 ,
147+ gaussPoint2
148+ ) ;
149+ let basisFunction = basisFunctionsAndDerivatives . basisFunction ;
150+ let basisFunctionDerivKsi = basisFunctionsAndDerivatives . basisFunctionDerivKsi ;
151+ let basisFunctionDerivEta = basisFunctionsAndDerivatives . basisFunctionDerivEta ;
152+ let xCoordinates = 0 ;
153+ let ksiDerivX = 0 ;
154+ let etaDerivY = 0 ;
155+ const numNodes = this . nop [ elementIndex ] . length ;
156+ for ( let nodeIndex = 0 ; nodeIndex < numNodes ; nodeIndex ++ ) {
157+ xCoordinates +=
158+ nodesXCoordinates [ this . nop [ elementIndex ] [ nodeIndex ] - 1 ] * basisFunction [ nodeIndex ] ;
159+ ksiDerivX +=
160+ nodesXCoordinates [ this . nop [ elementIndex ] [ nodeIndex ] - 1 ] *
161+ basisFunctionDerivKsi [ nodeIndex ] ;
162+ etaDerivY +=
163+ nodesYCoordinates [ this . nop [ elementIndex ] [ nodeIndex ] - 1 ] *
164+ basisFunctionDerivEta [ nodeIndex ] ;
165+ }
166+ for (
167+ let localNodeIndex = firstNodeIndex ;
168+ localNodeIndex < lastNodeIndex ;
169+ localNodeIndex += nodeIncrement
170+ ) {
171+ let globalNodeIndex = this . nop [ elementIndex ] [ localNodeIndex ] - 1 ;
172+ if ( side === 0 || side === 2 ) {
173+ // Horizontal boundaries of the domain (assuming a rectangular domain)
174+ residualVector [ globalNodeIndex ] +=
175+ - gaussWeights [ 0 ] *
176+ ksiDerivX *
177+ basisFunction [ localNodeIndex ] *
178+ convectionCoeff *
179+ extTemp ;
180+ for (
181+ let localNodeIndex2 = firstNodeIndex ;
182+ localNodeIndex2 < lastNodeIndex ;
183+ localNodeIndex2 += nodeIncrement
184+ ) {
185+ let globalNodeIndex2 = this . nop [ elementIndex ] [ localNodeIndex2 ] - 1 ;
186+ jacobianMatrix [ globalNodeIndex ] [ globalNodeIndex2 ] +=
187+ - gaussWeights [ 0 ] *
188+ ksiDerivX *
189+ basisFunction [ localNodeIndex ] *
190+ basisFunction [ localNodeIndex2 ] *
191+ convectionCoeff ;
192+ }
193+ } else if ( side === 1 || side === 3 ) {
194+ // Vertical boundaries of the domain (assuming a rectangular domain)
195+ residualVector [ globalNodeIndex ] +=
196+ - gaussWeights [ 0 ] *
197+ etaDerivY *
198+ basisFunction [ localNodeIndex ] *
199+ convectionCoeff *
200+ extTemp ;
201+ for (
202+ let localNodeIndex2 = firstNodeIndex ;
203+ localNodeIndex2 < lastNodeIndex ;
204+ localNodeIndex2 += nodeIncrement
205+ ) {
206+ let globalNodeIndex2 = this . nop [ elementIndex ] [ localNodeIndex2 ] - 1 ;
207+ jacobianMatrix [ globalNodeIndex ] [ globalNodeIndex2 ] +=
208+ - gaussWeights [ 0 ] *
209+ etaDerivY *
210+ basisFunction [ localNodeIndex ] *
211+ basisFunction [ localNodeIndex2 ] *
212+ convectionCoeff ;
213+ }
214+ }
215+ }
100216 } else if ( this . elementOrder === "quadratic" ) {
101217 for ( let gaussPointIndex = 0 ; gaussPointIndex < 3 ; gaussPointIndex ++ ) {
102218 let gaussPoint1 , gaussPoint2 , firstNodeIndex , lastNodeIndex , nodeIncrement ;
@@ -139,7 +255,8 @@ export class ThermalBoundaryConditions {
139255 let xCoordinates = 0 ;
140256 let ksiDerivX = 0 ;
141257 let etaDerivY = 0 ;
142- for ( let nodeIndex = 0 ; nodeIndex < 9 ; nodeIndex ++ ) {
258+ const numNodes = this . nop [ elementIndex ] . length ;
259+ for ( let nodeIndex = 0 ; nodeIndex < numNodes ; nodeIndex ++ ) {
143260 xCoordinates +=
144261 nodesXCoordinates [ this . nop [ elementIndex ] [ nodeIndex ] - 1 ] * basisFunction [ nodeIndex ] ;
145262 ksiDerivX +=
0 commit comments