Once we get all the parameters calculated and the matrix loaded, then we can finally solve the thing. Here’s how the code looks:
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 |
function loadMatrix(index,fceTemp) { var i; var j; var k; var tsa; var tfa; var ts3; var tf4; tfa = arrayTfce[index] + Tabs; tf4 = tfa * tfa * tfa * tfa; for (k = 0; k <= sizeNodes; k++) { ZSRCE[k] = 0.0; ZAUX[k] = 0.0; for (i = 0; i <= sizeBandWidth; i++) ZMATX[k][i] = 0.0; } k = 0; for (j = 0; j < nodesHigh; j++) for (i = 0; i < nodesWide; i++) { k++; //nodal energy ZMATX[k][1] += arrayRho[i][j] * arrayCp[i][j]; ZSRCE[k] += arrayRho[i][j] * arrayCp[i][j] * arrayT[i][j]; //conduction to right if (arrayKright[i][j] > 0.0) { ZMATX[k][1] += arrayKright[i][j]; ZMATX[k][2] -= arrayKright[i][j]; ZMATX[k+1][1] += arrayKright[i][j]; } //conduction up if (arrayKup[i][j] > 0.0) { ZMATX[k][1] += arrayKup[i][j]; ZMATX[k][1+nodesWide] -= arrayKup[i][j]; ZMATX[k+nodesWide][1] += arrayKup[i][j]; } //radiation if (arrayRad[i][j] > 0.0) { tsa = arrayT[i][j]; ts3 = tsa * tsa * tsa; ZMATX[k][1] += arrayRad[i][j] * ts3; ZSRCE[k] += arrayRad[i][j] * tf4; } } } //loadMatrix function oneCalc(index) { calcParams(index); for (var i=0; i<nodesWide; i++) { for (var j=0; j<nodesHigh; j++) { arrayT[i][j] = arrayTpieces[index][i][j] + Tabs; } } loadMatrix(index); calcMatrix(); k = 0; for (var j=0; j<nodesHigh; j++) { for (var i=0; i<nodesWide; i++) { k++; arrayTpieces[index][i][j] = ZSRCE[k] - Tabs; } } } //oneCalc var loopCount = 0; var s="Temperatures after "+loopCount+" time steps<br />"; for (var j=nodesHigh-1; j>=0; j--) { for (var i=0; i<=nodesWide-1; i++) { s += arrayTpieces[0][i][j].toPrecision(7) + ", "; } s += "<br />"; } document.getElementById("stuff").innerHTML = s; function runOneTime() { if (loopCount == 0) { for (var i=0; i< pieceCount; i++) { for (var j=0; j<nodesWide; j++) { for (var k=0; k<nodesHigh; k++) { arrayTpieces[i][j][k] = 70.0; } } } } else { oneCalc(0); } var s="Temperatures after "+loopCount+" time steps<br />"; for (var j=nodesHigh-1; j>=0; j--) { for (var i=0; i<=nodesWide-1; i++) { s += arrayTpieces[0][i][j].toPrecision(7) + ", "; } s += "<br />"; } document.getElementById("stuff").innerHTML = s; loopCount++; if (loopCount > 20) { loopCount = 0; } } var intervalID = setInterval(runOneTime, 1000); |
We start with a 7×7-node workpiece where the interior nodes are two inches wide and two inches high. The surface nodes are half as wide or high, meaning the cross-section area of each node is four square inches in the interior of the piece, one square inch on the four corners of the piece, and two square inches for the other surface nodes. The length of the surface exposed to the furnace is zero for interior nodes, zero for the bottom middle nodes (the piece is sitting on a solid hearth through which it is assumed no heat transfer occurs), one inch for the bottom nodes on each side, and two inches for the remaining surface nodes (the two corner nodes include one inch across the top surface and one inch down the side surface). The center of each interior node is assumed to be in the middle of the node, one inch from each side. The “center” of each surface node is assumed to be on the middle of the outer surface of the node. The “center” of the outer corner nodes is assumed to be on the outermost corner. This way there are always exactly two inches between the centers of each node, so the heat conduction terms make sense and the radiation can “see” the “center” of each node.
I have done it this way from the beginning and found that the readings of surface temperature on pieces which were actually heated in furnace matched the predicted values with impressive accuracy.
If all nodes are initialized to 70 degrees F, then exposed to a furnace at 1800 degrees F for successive 40-second time steps, we get the following results: