THE-15‘s
\ cyLml.‘-.;Q‘;‘ _, ... ..._ ... A... Iva-nulls-..
Ca \ ‘ i ‘1 I
I I
l “ 9 ~ -.'= a w
£15»; .sfl‘fit Yul ‘k i
ff?" as a A «A A ‘F;‘ f (7‘
i ‘- ‘-’~.-..a~ .5 $5.12." - V -w.‘ ‘
c 1-? -7 “ Q5.” "x“
‘1'“ 19b“ .3 we “#3,;
L '

 

 

This is to certify that the
thesis entitled
(omputer Simulation of Two-Dimensional

Transient Temperature Field in

Cryomicrosc0pe Conduction Stage
presented by

Siu-Ming Tu
has been accepted towards fulfillment

of the requirements for

Master (mgeeh, Mechanical Engineering

(Warm

a r professor

 

Date JUlY 19: 1983

0.7639 MS U is an Affirmative Action/Egan! Opportunity Institution

 

 

 

 

 

 

 

 

)V1ESI.] RETURNING MATERIAL§z
Place in book drop to
LIBRARIES remove this checkout from
.—:_—_ your record. FINES will 1‘
——- be charged if book is ,
returned after the date
stamped below. 5

 

2633M USE CL?

-' .

 

 

.fi I m" an!” — . , , N .-
7 . .: {V151, «sh gnu e-r- ‘ I 3 AW; ‘5'
- ' i u . 3

"Ht-rams diam

r*”s
t“ 3
FA.

5

t

I

a

F

G

I

(c

 

 

 

 

0/3 “$129352

COMPUTER SIIIJLATION OF ITO-DIMENSIONAL
TRANSIENT TEMPERATURE FIELD IN
CRYOMICROSOOPE CONDUCTION STIGE

BY

Sin-Ming Tu

A THESIS

Submitted to
lichigan State University
in partial fulfillment of the requirements
for the degree of

MASTER OF SCIENCE

Department of Mechanical Engineering

1983

ABSTRACT
COMPUTER SIMULATION OF TWO-DIMENSIONAL
TRANSIENT TEMPERATURE FIELD IN
CRYOMICROSCOPE CONDUCTION STAGE

BY
Sin-ling Tu

Computer controlled heat transfer systems have been developed for
microscOpes in order to allow visualization of the thermal response of bio-
logical cells to well-defined temperature changes.

Tbmperature measurement errors may arise in such heat transfer systems
due to the presence of thermal gradients within the sample and the fact
that only a single temperature sensor is used for the feedback control
100p. Since accurate temperature measurements are required over a wide
range of temperatures and a wide band width of cooling and warming rates. a
computer model has been deve10ped to simulate the response of such a sys-
tem. The model analyzed is a rectilinear composite material system with an
embedded plane heater and arbitrary thermal boundary conditions. The
two-dimensional transient response of this system is simulated using
interactive computer graphics.

The model was valid judging by comparison with experimental data. A
parameter study has defined the relative importance of various possible
design changes which would minimize temperature measuring errors and main-

tain a wide Operating bandwidth.

ACKNOWLEDGEMENTS

I would like to express my gratefulness to my advisor. Dr. J. J.
McGrath. who supported me and encouraged me during the course of this work.
I also would like to thank Dr. L. J. Segerlind (Agricultural Engineering
Department) for his permision of me to use his steady-state finite element
computer program. I also appreciate the assistance from Ray Thompson and
Guy Allen. They are consultants of the Case Center of Michigan State
University. and their assistance in the computer work of this research is
helpful.

Special appreication is due to James E. Otis. whose significant con-
tribution in the investigation of parameter variation solidifies this

research.

ii

TABLE OF CONTENTS

LIST OF TABLES......................................................
LIST OF FIGURES.....................................................
NOMENCLATURE........................................................
SUBSCRIPT NOMENCLATURE..............................................
GREEK NOMENCLATURE..................................................
CHAPTER

1. INTRODUCTION..................................................

2. EXPERIENTAL APPROACH FOR DETERMINING THE TEMPERATURE
DISTRIBUTION WITHIN THE CRYOMICROSCOE CONDUCTION STAGE......

2.1 Simplified Cryomicroscope Conduction Stage..............
2.2 Experimental Configuration and Procedure................
2.2.1 Purpose of Experiment.............................
2.2.2 Experimental Configuration........................
2.2.3 Experimental Procedure............................
2.3 Experimental Results....................................

3 . NUMERICAL APPROACH FOR PREDICI’ING THE TEMPERATURE DISTRIBUTION
WIND] THE CRYOMICROSCDPE CONDUCTION STAGE. . . . . . . . . . . . . . . . . .

3.1 Computer Model of Cryonicroscope Conduction Stage.......

3.1.1 Assumptions.......................................

3 .1 .2 Graphical Representation of Computer Model . . . . . . . .

3.1.3 Parameter Values..................................

3.2 Computer Program........................................
3.2.1 Program for Steady-State Temperature

Distribution....................................

3.2.2 Program for Transient Temperature Distribution....

4. COMPARISON AND DISCUSSION OF THEtEXPERIMENTAL AND NUMERICAL

mall“.O...I...OO...OOOOOOOIOOOOOOOOOOOOOOOOO...00.0.0000...

4.1 omparison of the Experimental and Numerical results....
42 i‘cnssioneeeeeeeeeeeeeeeeeeeeeeeeeeeeeeaeeeeeeeeeeeeeee

C
D
4.2.1 Discussion of Results Comparison..................
4
4

.202 Error Discn‘sion.O...OOOCOOOOOOOOOOIOOOOOO00......
.2.3 Discussion of the Graphical Output................

5. EFFECTS OF PARAMETER VARIATIONS ON THE THERMAL RESPONSE OF A
iii

V

vi

xii

xiii

12
12
16
19
20

22

22
22
24
28
33

33
37

58

58
85
85
87
90

“YOMImosmPE mNDUa‘ION STMEOOOOOOOOOOOOCOOOO...0.0.0....
5.1 Results of Parameter Variations.........................
5.2 Discussion and Conclusion...............................
5.3 Suggestions for Future Study............................
APPENDICES
A: C.1cu1‘tion Of P‘ruotct V‘lu‘...00....OOOOOOOOOOOOOOOOOOOOOO...
B: Derivation of Forward-Difference Equation for Interior Node......

C: Flowcharts and Listing of Cbmputer Program.......................

RMMQS.eeeeeeeeeeeeoeeeeeeeeeeeoaeoeeeeeeeeeeoeeoeeeeeeoeeeeoeee

iv

102
102

118
131

134

148

153

213

LIST OF TABLES

Experimental Conditions..........................................
Boundary Conditions..............................................
Thermal Properties and Heat Generation...........................
Heat Transfer Coefficients and Ambient Temperature...............
Quantitative Evaluation of Result Comparison.....................
Results of Investigation of Parameter Variations.................
Thermal PrOperties of Materials for Metal Plate..................

no SM‘ry 0f P‘rmeter vnri‘tionSOOOOOOOOOOOOOOOOOOOOOOOOOOOOOO

26

31

84

106

108

119

10.
11.
12.
13.
14.
15.
16.

17.

18.

19.

20.

21.

LIST OF FIGURES

crYMicroscope conduction st‘BOOOCOOOOOOOOOOOCCOOOOOOOOO0.0.0....
Simplified Cryomicroscope Conduction Stage (To Scale)............

Detailed Scheme of CryomicroscOpe Conduction Stage
(Not To sc‘l.)0..0..0.0.0.0.0000...OOOOOOOOIOOOOOOOOOOO0.0.0...

Experimental Configuration.......................................
A Typical Experimental Result of the Temperature History of
T1 on the Simplified Cryomicroscope Conduction Stage
(Beating Only).................................................
Graphic Representation of Computer Model.........................
Finite Element Grid for the Computer Model.......................
Result Comparison of Steady-State Temperature Distribution.......
Interior Node (with or without an Interface at the Node).........
Boundary Nodes..................................................
Corner Nodes....................................................
Inner Corner Nodes..............................................
Forward-Difference Grid for the Computer Model..................
Numerical Result of Nodal Temperature Distribution. . . . . . . . . . . . . .
Two-Dimensional Isotherm Plots (Case 1).........................

Normalized Three-Dimensional Temperature Distribution (Case 3)..

Superposition of Normalized Three-Dimensional Temperature
DIstribution at Different Time Instants (Case 3)..............

Tbmperature History of T; (Case 1)..........................

Tbmperature Difference between T1 and T: as a Function
OfTime (C‘se1)O0.0.000...OIOCOOOOOOOOOO0.00...0.00.00.00.00.

T1(t) OfC‘so100......0....OOOOOOOOOOOO0.000000000000000000

T3(t) OfC‘se10.0.00000000000000OO...OOOOOOOOOOOOOOOOOOOOOO
vi

17

21

25

34

36

40

41

42

43

49

51

52

53

54

56

57

59

60

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.

T,(t)
I;(t)
13m
T1(t)
T,(t)
T,(t)
am
nu)
ma
T,(t)
T,(t)
13m
T,(t)
mt)
T,(t)
13m
T;(t)
T,(t)
T;(t)
T,(t)
13m
mt)

mo

of

of

of

of

of

of

of

of

of

of

of

of

of

of

of

of

of

of

of

of

of

of

of

Case 1.............................................
Case 1.............................................
Case 1.............................................
Case 2.............................................
Case 2.............................................
Case 2.............................................
Case 2.............................................
Case 2.............................................
Case 3.............................................
Case 3.............................................
Case 3.............................................
Case 3.............................................
Case 3.............................................
Case 4.............................................
Case 4.............................................
Case 4.............................................
Case 4.............................................
Case 4.............................................
Case 5.............................................
Case 5.............................................
Case 5.............................................
Case 5.............................................

c.8° SOOOOOOOIOOOIOOOOOIOOIOOICOOOOOOOOOOOOOOOOOOOO

61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82

83

Transient Response of Isotherms (Cases 1 and 2)................. 92

Transient Response of 3-D Temperature Distribution

(Cases 1 and 2).

vii.

93

47.
48.
49.

50.

51.
52.
53.
54.
55.
56.
57.
58.
59.

60.

61.

62.

63.

64.

65.
66.
67.

68.

Transient
Transient
Transient

Transient
(Cases 4

Parabolic
T}(t) for
Effect of
Effect of
Effect of
Effect of
Effect of
Effect of

Effect of

Response of Isotherms (Case 3)........................
Response of 3-D Tbmperature Distribution (Case 3).....
Response of Isotherms (Cases 4 and 5).................

Response of 3-D Temperature Distribution

and 5)...............................................
Temperature Distribution within the Viewing Hole......
Investigating of Parameter Variation..............

Changing Width of Viewing Hole........................
Changing Thickness of Heating Glass...................
Changing Magnitude of Heat Generation.................
Changing Width of Heat Generation Region..............
Changing Thickness of Copper Plate....................
Changing Material of Metal P1ate......................

Changing Position of Coolant Channel..................

Thermal Resistance Diagram for a Composite Material with
Heat Generation.00......0......0.0.000...OOIOOOOOOOOOOOOOOOOOO

Comparison with Khompis's Result of Changing the Width of Heat
Generation Region (Khompis's Definition for Normalization)....

Comparison.with Khompis's Result of Changing the Thickness of
Glass Heater (Khompis's Definition for Normalization).........

Comparison with Khompis's Result of Changing the Width of Heat
Genet‘tion Region.OOOCOOOOOOIOOOOOOO0.00.00.00.00...00.0.00...

Comparison with Khompis's Result of Changing the Thickness of

61‘83 Ho‘tetOOOOI0..OO...00....0.0.0.0000...0.00.0000...00....

Descriptive Scheme of the Recommendation........................

Extrapolation for Thermal Conductivity of Kapton Film...........

P.th 0f coal‘nt FIWOOOOOOO....GOOOOOOOOOIOOOOOOOOOIO0.0.00.0...

Interior Nodal System with Energy Terms Indicated...............

viii

95

96

98

99

100

104

111

112

113

114

115

116

117

122

126

127

128

129

132

136

143

149

69.
70.
71.
72.
73.
74.
75.
76.

77.

Divisions

Flowchart

Flowchart

Flowchart

Flowchart

Flowchart

Flowchart

Flowchart

for I'othon Plottin‘...0.0...OOOOOOOOOOOOOOOOOOOOO...

of

of

of

of

of

of

of

“in PrO‘rm 'BET'OOOOOOOOOOOOOOOOOOOOOOOOOOOOOOO

Subroutine

Subroutine

Subroutine

Subroutine

Subroutine

Subroutine

'INOUT'.................................
'GRAPH3D'...............................
'GRAPHZD'...............................
'ISOTHERM'..............................
'SIDEB'.................................

'm'OOOOOOOOOOOOOOCOOOOOOOOOOOOIOOOOO

hl‘tion “on: Prosrm Unit‘OOOOOOOOOOOOOOOOOOOOOCOCOOOOOOOOOOOO

ix

155

157

160

162

163

165

169

170

171

Nu

Pr

Ra

Re

Rt

RY

NOMENCLATURE

Area
Specific Heat
Normalized Heat Generation Per Unit Area

Heat Generation within Nodal System. or Gravitational
Acceleration

Grashof Number

Heat Transfer Coefficient

Thermal Conductivity

Length

Mass Flowrate

Nusselt Number

Pressure. or Power Input

Prandlt Number

Heat Flux

Electrical Resistance

Rayleigh Number

Reynold Number

Dimension Ratio. Lx/I.y

Thermal Conductivity Ratio. K3]!1

Time Constant Ratio. too/(L!)a

Ratio of Grid Size in XrDirection. szle1
Ratio of Grid Size in Y-Direction. AYa/AY1
Tempe rature

Time

Electrical Voltage across the Glass Heater. or Flow Velocity.
or Volume

X—Coordinate
Y-Coordinate

Elevation of Coolant Flow

xi

atm

bo

SUBSCRIPT NOMENCLATURE

Silver Epoxy

Atmosphere

Bottom

Bottom Surface of the Viewing Hole

Coolant

Glue Epoxy

Kapton Film. or Quantity Evaluated at Film Temperature
Glass

Nodal-Numbering in XrDirection. or Imaginary Material
Nodal-Numbering in Y-Direction

Left

Quality

Right

Top

Wall

X-Direction

Y‘Direction

Ambient Characteristics

SUPERSCRIPT NOMENCLATURE

Time Interval

xii

GREEK NOMENCLATURE

Thermal Diffusivity
Difference

Non-Dimensional Temperature
Density

Dynamic Viscosity

xiii

CHAPTER 1

INTRODUCTION

Many important medical problems are related to the heat transfer
characteristics of cells and tissues. Examples of such problems are:
the selective cellular destruction by freezing to minimize blood loss.
the killing of tumors. and alleviating recovery pain by the temporary
destruction of nerves. Generally speaking. cells and tissues are quite
sensitive to temperature. temperature change. and the rate of tempera-
ture change. For example; the survival rate of CHO cells during the
heat-induced killing process will change from 10% to 0.1i when the heat-
ing temperature °h338¢8 from 420C to 42.5'C (reference [1]). Hence
temperature differences as small as one half degree Celsius significant-
ly affect the survival rate of CHO cells. Therefore. it is worthwhile
to study and quantify the response of cells and tissues to various tem-

peratures and rates of temperature change.

2

A temperature controllable device is needed to study the heat
transfer characteristics of cells and tissues. Several types of tem-
perature controllable device have been deveIOped and used in the
laboratory. One of these devices is the cryomicroscope conduction stage
which is used compatibly with a microscope to observe the behavior of
biological samples under different temperature environment and different
rate of temperature change. A cryomicroscope conduction stage was
designed by Dr. J. J. McGrath as a modification of the
Diller-Cravalho cryomicroscOpe stage (reference [2]). The purpose of
this modification was to allow fluorescence optics to be employed as a

mean to detecting cell viability on the microscope.

Figure 1 shows the cryomicroscope conduction stage used in the
Bioengineering Transport Processes Laboratory within the Mechanical
Engineering Department of Michigan State University. Basically. the
system is composed of microscOpe cover glasses. a circular quartz
heater. and a copper plate. A coolant channel is provided inside the
capper plate. The coolant is introduced into the coolant channel and
circulates around the circular quartz heater to provide a cooling
effect. On the other hand. a thin layer of tin oxide coated on one side
of the circular quartz disc such that it becomes a transparent electric
resistor. When electric current is applied heat is generated in the tin
oxide coating which provides a heating effect. Thus. the system is

cooled and heated simultaneously.

The existing system manipulates the temperature and rate of tem-

perature change around the viewing hole by proper control of the heat

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

. \\
. / / \ \ .
“D \ QL/_il \:_X”mu E
"f1. (EV/776D" {33
©——— \5: " i ’ -——®

 

 

 

(9

1. COVER GLAss 2. CIRCULAR QUARTZ HEATER 3, pLASTIc CLAMP
4. NYLON BOLT 5. MYLAR TAPE 5. METAL CONDUCTOR
7. VIEWING HOLE 8. COOLANT CHANNEL 9. COPPER PLATE

lO. ELECTRIC POWER WIRE SLOT 11. COOLANT INLET AND OUTLET

FIGURE 1. CRIOMICROSCOPE CONDUCTION STAGE

4

generation in the quartz disc. keeping the coolant effect constant.
Thus. biological samples can be placed on the viewing hole so that the
behavior of such samples can be observed under the microscope in the

real time as the temperature of the sample is changed.

Two pieces of cover glass are bonded on the tOp of the circular
quartz disc so that a temperature sensor can be embedded between the two
cover glasses without contacting either the biological samples or the

heating surface.

Plastic clamps keep the circular quartz heater in contact with the
metal conductor which conducts electric current into the tin oxide coat-
ing. In the present case. the tin oxide is coated on the lower surface

of the circular quartz disc.

A thin layer of silver epoxy also has been painted onto the two
Opposite sides of the lower surface of the circular quartz disc so that
the electric current passing through the tin oxide thin film will be
more uniform. Mylar tape is placed on top of the cOpper plate to pre-

vent short circuiting of the thin film through the copper plate.

Generally, an ideal cryomicroscope conduction stage should possess
an uniform temperature distribution within the viewing hole where sam-
ples are observed. because only one temperature sensor is used to
measure the temperature of the biological sample. Thus. the temperature
measurement error will be reduced. if the temperature distribution with-

in the viewing hole is uniform.

5

An ideal cryomicroscope conduction stage also should be capable of
Operating over a wide range of rates of temperature change. Optimal
freezing rates for various cell types range from approximately 0.1
oC/min to 10.000°C/min which suggests that such rates will be required
if these cell types are to be examined with such a system. These two
characteristics become the criteria of the design of cryomicroscope con-
duction stage. That is. an optimal design of the cryomicroscope
conduction stage is required . which can minimize the measurement error
due to treating the temperature of the temperature sensor as ‘the tem-
perature of the biological samples and can be operated under a wide

range of thermal response.

Heat transfer analysis of the cryomicroscope conduction stage has
been done by Vorapot Khompis (reference [3]). Khompis has built a com-
puter model for the cryomicroscOpe conduction stage and used the
forward-difference method to solve the temperature distribution within
the system. The present research can be considered as a continuing
study of his work in more detail of the computer model and more work on

the study of the parameter variations.

Basically. the purpose of this research is to simulate the dynamic
response of the temperature distribution within the cryomicroscOpe con-
duction stage using a computer heat transfer model. Using this model.
the effects of various design parameter variations on the uniformity of
the temperature distribution within the viewing hole and the rate of
temperature change at the position where the biological samples are

rested are investigated. These computer simulations are then used to

6

recommend an improved cryomicroscope conduction stage design. Prior to
use of the computer model to study the effects of parameter variation
the model was verified by comparison with experimental temperature his-
tories obtained from the cryomicroscope system in a number of

experimental configurations.

CHAPTER 2
EXPERIMENTAL APPROACH FOR DETERMINING THE TEMPERATURE DISTRIBUTION

WITHIN THE CRYOMICROSCOPE CONDUCTION STIGE

2.1 SIMPLIFIED CRYOMICROSCOPE CONDUCTION STAGE

In the interest of obtaining a more accurate computer simulation.
the original experimental cryomicroscope system has been modified to
more closely match the necessary simplifications made in the computer
model. Figure 2 shows the simplified cryomicroscope conduction stage.
Figure 3 shows the detailed schematic cross-sectional area of the sim-

plified cryomicroscOpe conduction stage.

Several differences between the original experimental
cryomicroscope system (Figure 1) and the simplified system (Figure 2)

are stated as followings:

 

 

 

 

 

P
l}.

/ 9

1L....
n"---
H

II

/

 

 

 

 

 

 

 

 

 

6
”G.--

 

 

 

 

 

 

 

 

 

 

 

 

 

I:
r
‘ljg

1. COVER GLASSES 2. METAL CONDUCTOR 3. PLASTIC CLAMP

 

U/@
a

h"-—(
Hulk-4

f 8
a “SI—(7)

l
I
,—1
l

 

 

HU—A

 

 

 

4. RECTANGULAR HEATING GLASS 5. VIEWING SLOT 6. COPPER PLATE
7. COOLANT CHANNEL 8. ELECTRIC POWER WIRE SLOT
9. BOLT FASTENER

FIGURE 2. SIEPLIFIED CRYOKICROSCOPE CONDUCTION STAGE (TO SCALE)

Amq<om 09 902v dflm< Adonaommlmmomo.mmoomOGOHzowmu ho.mxmmom DNAHdfimn .n RQDOHh

A ozHa<oo HQHNO DHAA<HHZ v :AHm m>Hsoonzoo NAA<OHmBofiAm or; mz<qo UHam<Am on—

aoqm mmHB mmzom uHmaomAm .N— Amzz<mo Bz<qooo .p. moz<m¢mgo ao<szoo .o— scam oszmH> .m

quh zoam<s .Q flags on ma<qm mmmmou .m mm<qo ozHB<flm on
axomm mH>AHm .: moeu=nzoo A<Bmz .n mm<do mm>oo .N HDAO on

a. g .3 a

 

 

 

 

 

 

 

 

 

 

 

 

7.;
a,

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

lt/@

x

 

 

 

 

10

A. Rectangular Glass Heater: One of the basic assumptions of this
research is that the temperature distribution is uniform along the
Z-direction (see Figure 2). that is. only the temperature distribution
in the XFY plane is of concern. Tb fulfill this assumption. the circu-
lar quartz heater in the original system has been replaced by a
rectangular glass heater in the new system.

In the original system. 'square' cover glasses are placed on the
'circular' quartz heater which. in turn. is placed on the 'rectangular'
capper plate. This situation creates a three-dimensional problem
because of its combined geometric configuration. Tb simplify the prob-
lem. the basic configuration of the cryomicroscope conduction stage has
been modified to be rectilinear so that it can be treated as a

two-dimensional heat transfer problem.

B. Heating surface: The surface coated with the metallic oxide of
the rectilinear glass heater in the new system is installed upward.
while in the original system the heating surface is faced downward.

A thin layer of Mylar tape which had been placed on the copper plate to
prevent the copper plate from electric contact with the metallic oxide
coating can be removed as a result of moving the heating surface.
Further. two pieces of metal conductors which had been inserted between
the Mylar tape and the circular quartz heater to serve as electric ter-
minals for conducting electric current into the metallic oxide coating
can also be moved to the upper surface of the glass heater. This rear-
rangement greatly reduces the possibility of short circuits especially
when ice is growing on the cryomicroscope conduction stage during exper-

iment.

ll

In the original system. the clearance between the circular quartz heater
and the Mylar tape due to the thickness of the metal conductor created a
barrier to heat transfer between them. Since in the new system the
metal conductors are placed on the upper surface of glass heater. there
will be no significant clearance between the glass heater and the capper
plate. There will be a thermal contact resistance created at this junc-

tion however.

C. Position of Coolant Channel: The coolant channel has been
reconstructed so that only the two opposite sides of the copper plate
are involved in the heat exchange with the coolant channel (see Figure
2).

In the original system. the coolant channel was built inside the copper
plate and surrounded the circular quartz heater. This configuration

increases the expense of the system.

D. Viewing Slot: The circular viewing hole of the original system
has been changed to a rectangular slot to fulfill the assumption of a
two-dimensional heat transfer problem where the temperature distribution

along the Z-direction is considered as uniform.

E. Kapton Film: Since both the glass heater and the copper plate
are essentially rigid. the contact area between them becomes nonhomo-
geneous and unpredictable . which reduces the reproducibility ~of the
experimental results.

Also. nonhomogeneous contact will affect the temperature distribution on

the top surface of cover glass because the thickness of the glass heater

12

together with the cover glasses is relatively small compared to the
other dimensions of the system (see Figure 2).

To improve this situation. a thin. flexible Kapton film is inserted
between the glass heater and the copper plate to reduce the thermal con-
tact resistance by increasing the contact area and diminishing the
indefiniteness of the contact effect. The Kapton film is a product of
Du Pont Company, and it has the ability of maintaining its properties

over wide temperature range.

F. Dimension Change: The dimensions of the system have been
changed relative to the original design. and this is subjected to an
idea that setting a standard size for different microscope stages such
as concentration diffusion chamber and cryomicroscope conduction stage
so that they all can match with a standard frame which is used compati-

bly with microscope to hold the stage at the right position.

2.2 EXPERIMENTAL CONFIGURATION AND PROCEDURE

2.2.1 PURPOSE OF EXPERIMENT

To simulate the dynamic behavior of the temperature distribution of
the cryomicroscOpe conduction stage. several thermocouples were embedded
at certain important locations to measure the temperature histories at
these locations. This discrete information then was used to represent
the dynamic behavior of the temperature distribution of the cryomicro-

scope conduction stage. Since in the simplified system the temperature

13
distribution along the Z-directicn is neglected. these thermocouples
were distributed within the cross-sectional area of the system. that is.

in the X-Y plane.

Generally, it is impractical to obtain experimental or
computer-simulated temperature histories for every point of the system.
so it is reasonable to get temperature histories at certain important
locations which should be carefully selected. The selection should be
such that if the temperatures at these selected points are matched by
those of the corresponding locations from the numerical results. then
the temperature distribution of the cross-sectional area of the cryomi-
croscope conduction stage can be assumed to be described by the

temperature distribution from the numerical results.

Figure 3 indicates the locations selected for embedding thermocou-
ples for this system. They are T1, T3, 1., 1;, ‘nd Ts, T; is the
location where the temperature sensor is usually embedded to measure the
temperature of biological samples. In this research. no thermocouple
was placed at this location. The positions of T1 and T, are selected so
that the temperature difference between them. AT. can serve as an indi-
cation of the uniformity of the temperature distribution within the

viewing hole area.

The difference between T1 and T, gives an indication of the tem-
perature difference across the thickness of the system. Ti, T3, and T,
reveal essential information concerning the temperature distribution on

the top surface of the system. Since the thermal conductivity of copper

11+
is significantly greater than the other materials in the system. it is
assumed that the temperature distribution of the copper plate is quite
uniform. Thus. only one thermocouple is placed in it. The numerical
results show that there is about a 5°C difference across the copper
plate when a total temperature range of 200’C occurs within the entire

system during operation.

Furthermore. this amount of difference will only have a trivial
influence on the temperature distribution of the region of concern.
namely. the viewing hole. because they are far apart. If the tempera-
ture histories at these selected locations for the experimental results
and the numerical results match. then the computer model is assumed to
be valid.

Several cases are considered in order to assess the general validi-
ty of the computer model rather than assuming that validity is
established by a match of experimental and predicted results for a sin-

gle case. Table 1 lists the various cases considered.

The reason to test the computer model under these conditions is to
make sure that it is not only valid for simple heating or cooling condi-
tions (Cases 1. 2. and 3) but that it can also predict the results in
the cases of superpositioning of heating and cooling (Cases 4 and 5).
Due to the linearity of the governing equation (2.1) for two-dimensional
heat conduction. superposition is expected and can serve as a check on

both experimental and numerical results.

Militias
Heating
Cooling

Voltage Across
Glass Heater

Resistance of
Glass Heater

Ambient Temp.:
Tt(Top)
T§(Bottom)

Tc(Coolant)

uNo coolant flow

TABLE 1

EXPERI DENT AL mNDITIONS

Yes

No

7.54

51.

23

23

IH

No

8.47

51.7

23

23

Ease}. Cassi Que;

No

Yes

23

-100

~196

Yes
Yes

7.66

54.

23
-100

-196

Yes
Yes

8.61

54.

23
-100

-196

Dimension

Volt

Ohm

16

8’T 3’1“ 1 u
+ = —---—— (2.1)
ax’ ar’ e at

 

 

However. the actual system is not expected to behave exactly like a
linear system. because the thermal properties will generally be tempera-

ture dependent over such a wide temperature range.

2.2.2 EXPERIMENTAL CONFIGURATION

Figure 4 represents a schematic configuration of the experiment.
Coolant available from the Liquid Nitrogen Tank passes through a Control
Valve and is introduced into the coolant channel of the Cryomicroscope
Conduction Stage. It then flows out of the stage to the Plenum. and
discharges through the orifice in the wall of the Plenum. The flowrate
is calculated from the total pressure which is measured by a pressure

tap and U-Tnbe (see Appendix A).

As soon as the switch is turned on. electric current passes through
the metallic oxide coating on the glass heater which then generates heat
by its electrical resistive characteristics. The temperature histories
at selected points are recorded by a Strip Chart Recorder. Since only
one Strip Chart Recorder is used here. the anaIOg signals from different
thermocouples are displayed periodically on the Strip Chart Recorder

using a Multi-Switch which is operated manually.

The information on strip chart was then discretized by eye and

input to the PRIME 750 Computer which was used to convert these input

1?

BDABDO
UHmmmmw

onadethzoo AdBZHZHmmmxH .: mmeHh

mmnmoomm
Bm<mo mHmBm moBHgthBADE

 

 

 

 

 

mmset:

 

 

 

 

 

 

 

 

 

 

zmsmmm
$9328 1 .. .. u ..
emu mszm /

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

xz<9
zmoomBHz
QHDGHA
% .L
m>q<>
domazoo

 

mmam am
ozHaooma
mam mmHa
momaooozmmme
232mmm
moHaHmo
may
mmammmmm 15(s(>(i

 

 

 

mu<em ZOHHODQZOO
mmoomomomzowmo

 

wqmmam
mmaom
on

..T

 

 

moaHam -\\\v

 

 

18

data into a graphical output suitable for comparison with the computer

simulation results generated on the PRIME 750.

Another Strip Chart Recorder and thermocouple was used to monitor
the coolant temperature . because the dynamic temperature response of
the system is quite sensitive to the coolant temperature. The individu-
al components are described as followings:

A. Simplified Cryomicroscope Conduction Stage: The electrical resistor
represents the metallic oxide coating.

B. Multi-Switch: For performing periodic input to the Strip Chart
Recorder.

C. Strip Chart Recorder: For recording temperature histories at
selected locations of the system.

D. Liquid Nitrogen Tank: For supplying the coolant which is a
vapor-liquid mixture. The vapor comes from fast evaporation during the
flow process.

E. Control Valve: For controlling the flowrate of coolant flow.

F. Plenum: From which the total pressure of the coolant flow is meas-
ured so that the flowrate of coolant can be calculated (see Appendix A
for calculation).

G. U-Tube: Measures total pressure of the coolant flow.

R. Variable DC Power Supply: Supplies electric power for heat genera-
tion.

1. Switch: For creating a step input for heat generaton.

J. PRIME 750 Computer System: Available in the Case Center of the
Mechanical Engineering Department of Michigan State University. This is

used to convert the raw data into a comprehensive graphical output for

19
comparison with computer generated model predictions.
K. Thermocouple: For measuring temperature histories at selected loca-
tions of the system and ambient temperatures. In addition. a digital
multimeter is required to measure voltage and electric resistance of the
metallic oxide coating. and a thermometer is used to monitor the ambient
temperatures and coolant temperature in the Plenum which is required in

calculating the flowrate.

2.2.3 EXPERIMENTAL PROCEDURE

A. ann on Variable DC Power Supply and Strip Chart Recorder.

B. Take the reading from the digital multimeter for the resistance of
the metallic oxide coating. Set the voltage magnitude across the metal-
lic oxide coating for open-100p conditions so that once the Switch is
on. and the circuit is activated . there will be a step input of heat
generation.

C. Open the Control Valve and/or turn on Switch.

D. Start turning Multi-Switch manually for periodic input of different
thermocouples. until the system reaches its steady-state.

E. Take the reading from the thermometer for the ambient temperatures
above and below the system and coolant temperature inside the Plenum by
thermocouple measurement.

F. Thke the readings of the water level at both columns of the U-Tube
so that the total pressure and flowrate can be calculated.

G. The experiment is complete.

20

2.3 EXPERIMENTAL RESULTS

A typical experimental result in graphical form as it exists on the
PRIME 750 is shown in Figure 5. Complete information will be given in

Figures 20 to 44 overlapped with numerical results from the computer for

the convenience of comparison.

21

Aquo wzHemmmv mwmem mmoomomonowmo QmHthmsz mme 20
PB Lo meBme ambedmmmzme Mme ho Bubmmm adazmszmmxm Q¢0Hm>8 ¢ .m mmeHm

.umm.

-.~m_ ..mou ..<<. ...~_ -.mm

o.~n -.m< ..<m c.-

 

 

4».

\

 

 

 

 

.IIAf-‘lvslsl Ill

A“

1“

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

ms.» w>

L

mmah<mmazm»

...~
..m~
o.-m
a.mm
o..<
a.m<
a.am
u.mm
..-o
a.mm
a..~

.u wwwmwo.

CHAPTER 3
NUMERICAL APPROACH FOR PREHCI‘ING THE TEMPERATURE DISTRIBUTION

WITHIN THE.CRYOMICROSCOPE CONDUCTION STAGE

3.1 COMPUTER MODEL OF CRYOMICROSCOPE CONDUCTION STAGE

3 .1 .1 ASSUMPTIONS

There are several important assumptions that have been made in the
development of the computer model to be presented in this section.

These are outlined at this time.

The first assumption is that the temperature distribution in the
Z-direction (see Figure 2) can be neglected. that is. only the tempera-
ture distribution in the XrY plane is of concern. Therefore the present

problem is considered to be a two-dimensional heat transfer problem.

22

23
The second assumption is that if the numerical results at selected
points match that of the experimental results at corresponding loca-
tions. then the temperature histories at arbitrary positions of the
system can also be predicted by the numerical results. It is also
assumed that the behavior of AT. ( the temperature difference between T1
and T3). and the rate of temperature change of T; of the confirmed come

puter model will follow that of the actual system.

The third assumption is that the thermal conductivity and
diffusivity of each material will be constant with respect to tempera-
ture. This constant value will be selected as the average thermal

prOperty value of the operating temperature range.

Due to the fact that the cover glass is thin. the temperature
difference between T1 and T; is neglected. where T; is the location to
be put temperature sensor to measure the temperature of the biological
samples during the experiment (see Figure 3). From the
computer-predicted results of the computer model which is proved to be
valid in chapter 4 the temperature difference between T1 and T; is with-
in 1°C for most of the cases. A typical computer-predicted temperature
distribution can be seen in Figure 14 which contains the nodal tempera-

tures of the computer model.

It is also assumed that the boundary adjacent to the Plastic Clamp
used to clamp the glass heater to the copper plate (see Figure 3) is
adiabatic due to the poor conductivity of plastics. The heat transfer

due to radiation is assumed to be negligible also.

21+

The last assumption is that the metallic oxide coating in the cen-
tral portion of the glass heater between the silver epoxy generates heat

homogeneously (see Figure 3).

3.1.2 GRAPHIC REPRESENTATION OF COMPUTER MODEL

Figure 6 represents the physical system as it is modelled in the
computer program. It is a two-dimensional. composite material. mixed
boundary condition . heat transfer problem with regions of planar heat
generation. The basic configuration represents the cross-sectional area
of the cryomicroscope conduction stage which is perpendicular to the

Since this cross-sectional area is symmetric with respect to the
center line. it is possible to consider only one half of the system.
Use is made of an adiabatic boundary condition on the center line. The
boundary conditions for this computer model are listed in Table 2. and

the numbers used to indicate the boundary regions are shown in Figure 6.

Several points should be made at this time. The first point is
that in the actual system. the glass heater and the Kapton Film are
clamped together by a plastic clamp. so there is contact resistance

between them from the heat transfer point of view.

But in the computer model it is imagined that there exists a thin
layer of material between the Kapton Film and the lower surface of the

glass heater. and the thermal conductivity of this imaginary material is

5
2

:WJNOoO

:NOOoO
:MOOoO“
:—O0.0

:NNOoO

:FOOoOn
:QOOoOL
:FOOoO
:QOOoO

 

Ammo: mmezmzoo mo onB<Ezmmmmmmm onm<mo

:meOoO

T..m_.o+| .985 |+|+| ..mm.ol.+lll|. _.m.o IlliTliméll

@

n:

x

.w mmeHm

6

 

 

 

 

 

 

 

 

 

 

 

 

 

_

f

 

 

a _
.H
III by ‘ . x
u I .2392: 3235: m E EES
on: : ., mo mfifi <
I _ 1:0/ / Q. ‘ . ©
n© * r I \AWA\%\\\\\“ fl
- s . _
\ @ V 2 E; _ _
\\\ \\ . 0*; _
“ 282$sz magic mmiw W apes 0:25 _ S :25 _
SE \ .. .,.. m
\ ® @ \ _ 52508 _
T. \ . n _ @ . _
r v N 35; mat/8 _ _
_I r
\ :. main $38 It mmI
[AWW ~\ HT+
. EEO @ .
_

.I

Boundary

9--10--11

11--12

26

TIBLE 2

BOUNDARY CONDITIONS

andition

Free Convection on the Lower Surface of a Horizontal Plate
Free Convection on a Vertical Plate

Free Convection on the Lower Surface of the Viewing Hole
Adiabatic Condition due to Symmetry

Free Convection on the Upper Surface of a Horizontal Plate
Free Convection on a Vertical Plate

Free Convection on the Upper Surface of Horizontal Plate
Equivalent Heat Convection to account for the Effect of
Contact with Metal Conductor

Adiabatic Condition due to Assumption of Poor Thermal
Conductivity of Plastic Clamps

Forced Convection due to Coolant Flow

27
varied to take care of the effect due to the contact resistance between
the Kapton Film and the glass heater (compare Figure 3 with 6).
Otherwise. if the contact effect is neglected or it is assumed that the
Kapton Film and the glass heater are bonded together. then the computer
will predict that the temperature distribution between the Kapton Film

and the heating glass is quite smooth.

However. as expected in the actual situation. the experimental data
show that the magnitude of the temperature changes abruptly across the
gap between the glass heater and the Kapton Film (see Figure 16). This

temperature difference can be as large as 40°C (case 5 of Table 1).

The experimental results also show that inserting the Kapton Film
can minimize variations in the experimental results due to variable
thermal contact resistance between the glass heater and the Kapton Film.
This makes the experiment more reproducible. because the flexibility of
the film increases the contact area between the two rigid surfaces mak-
ing the heat transfer between them more homogeneous with respect to

space.

In the actual system (see Figure 3). two metal conductors are con-
nected to the silver epoxy to conduct electric current into the metallic
oxide coating . and experimental results show that the effect of these
metal conductors can not be neglected. During a heating process (case 1
in Table 1). the experimental result of T, (see Figure 3) which is near
the metal conductors is lower than the computer-predicted result if the

computer model does not take care of the effect of the metal conductors.

28

During a cooling process (case 3 in Table 1). the experimental result of
T, is higher than the computer-predicted result if the computer model
does not include the effect of metal conductors. After inspecting the
heat transfer characteristics of the system. it is only the metal con-
ductors that can affect this portion of the system in this way. So. a
particular heat transfer coefficient. H“ is used to take care of the

effect due to metal conductors. The value of H.8 is determined from the

process of matching the experimental and numerical results.

When coolant is introduced into the system. the system temperatures
are reduced. In particular. the ambient temperature beneath the metal
plate (see Figure 2) will also decrease significantly. because there is
only a small clearance between the lower surface of the metal plate and
the table upon which the system rests. This suggests the use of a dis-
tinct ambient temperature for describing the free convection from the
lower surface of the system. This temperature was measured experimen-
tally aud was found to remain constant. The value of the temperature

was -100°C and this value was input to the computer program.
Furthermore. the heat transfer coefficients for the lower surface

of the copper plate and that of the viewing hole will be distinguished

because of their different geometric configurations.

3.1.3 PARAMETER VALUES

In some respects. 'HDDELLING' is a process of data matching. Thus.

some of parameters in the computer model are fixed while others are

29

ill-defined. This latter group requires some best fit matching pro-
cedures. and this is done by adjusting certain parameters of the
computer model. However. not all of the parameters of the computer
model can be adjusted. because some of them are already documented. such
as the thermal properties of copper or the magnitude of the heat genera-
tion (which can be calculated reasonably and easily) (see Appendix A).
This type of parameter will be called 'nonradjustable parameters'.

because they can not be changed during the process of data matching.

Also there are parameters that are adjustable. such as the heat
transfer coefficients for heat convection on the boundary. because the
calculation for these coefficients is based on some rough assumptions.
It is therefore reasonable to have the privilege to modify these values
as long as it can bring the computer predicted results close to the
experimental results. This type of parameter will be called 'adjustable

parameters'.

In addition there are parameters that are difficult to evaluate and
have to be found from data matching by trial and error. such as the
thermal properties of the imaginary material which is used to take care
of the effect due to the contact resistance between the Kapton Film and
the glass heater. Since it is imagined to be present. it is impossible
to evaluate its thermal properties. This type of parameter will be

called 'undetermined parameters'.

A list of the three types of parameters is given below.

Non-adjustable Pgrameters:

30
thermal conductivity and diffusivity of copper. K0 and ac. thermal con-
ductivity and diffusivity of glass. KI and “8' thermal conductivity and
diffusivity of Kepton film. If and of. magnitude of heat generation. G.
Adjustable Pgrameterg:
heat transfer coefficient of upper surface. Ht. of lower surface. Hb. of

vertical surface of copper plate. H of vertical surface of cover

r'
glasses. 31' of coolant flow. Ho. of lower surface of viewing hole. nbo'
Ind H‘B' a particular coefficient mentioned before.

Undetermined Pgrgggters:

thermal conductivity and diffusivity of glue. Xe and “e' thermal conduc-

tivity and diffusivity of imaginary material. Xi and ai’

Since there are many adjustable and undetermined parameters. it
would be impossible to determine their proper values by trial and error
on the computer. The strategy used here is to use a computer program
for the steady-state heat transfer problem to match the steady-state
experimental results of Case 1 (see Table l). which does not involve the
cooling effect. This will eliminate the thermal diffusivities. which
include capacitive effects important in the dynamic temperature
response. and Re. which is the heat transfer coefficient for forced con-
vection due to coolant flow. After obtaining a set of satisfactory
values for the static parameters. then computer-matching of the dynamic

temperature response data was undertaken.

Tables 3 and 4 list these parameter values for each case. A
detailed description concerning the determination of each parameter

value will be presented in Appendix A.

Pgrgetgr

Operating
Temp. Range

31

TABLE 3

THERMAL PROPERTIES AND HEAT'GENERATION

Cgsg 1
(23.50)

403
1.436
0.158
0.389
0.0025
1.16E—4
8.27E—7
1.02E—7
0.2EP7
1.E#6

0.1724

.Dimensionless

£ass.1
(23.55)

403
1.453
0.158
0.389
0.0025
1.16E-4
8.27Er7
1.02E-7
0.2Ek7
1.E-6

0.203

£ase.i.

£a12.1.

£asa.£

(-196.23) (-196.26) (-196.26) “c

495
0.796
0.135
0.389
0.004
1.33EF4
8.76Ek7
0.94E-7
0.2E-7
1.E-6

0.0

495
0.95
0.136
0.389
0.004
1.33Ek4
8.6E—7
0.94E-7
0.2Er7
1.E-6

0.27

484.4
0.995
0.137
0.389
0.004
1.33E-4
8.6Ek7
0.94E-7
0.2EF7
1.E#6

0.3528

vlu-‘x
WIN-'K
vlu-‘x
IYH-°K
Ila-“x
u’ISec
H’lSec
I’ISec
M’lSec

l3/Sec

n
as
“so

Tt(Top)

52

TABLE 4

HEAT’TRANSFER COEFFICIENTS AND’AMBIENT TEMPERATURES

59

33.92

101.75

6.09

Tb(Bottom) 23

Tc(Coolant) “

..No coolant flow

32.79
11.57
97.4
59.3
32.79
98.37
5.94
23

23

14.75
95.2
54.3
578
22
14.75
23
-100

-196

52.3

583

21.5

14.75

23

-100

-196

55.1
676.8
21.5
14.07
23
-100

-196

Dimension

wvu‘o'x
llu’-°x
w7u’-°x
w7u‘-'x
vvu’-'x
v/u’-'x

Ilu’-°x

33

3.2 COMPUTER PROGRAM

3.2 .1 PRmRAM FOR THE STEADY-STATE TEMPERATURE DISTRIBUTION

This computer pragram is based on the finite element method using a
linear triangular element and a bilinear rectangular element. It can be
used to solve any steady-state. two-dimensional field problem. This
program was developed by Dr. Larry J. Segerlind (Agricultural
Engineering Department) for academic teaching in a finite element course
offered at Michigan State University. The output of this program
includes the steady-state temperature at each node and the temperature
gradient of each element. Since only the linear triangular element and
the bilinear rectangular element are used in the program. it is easy to
define the temperature gradient of an element by specifying the tempera-

ture gradient at the central location of this element.

Generally. this gradient indicates the orientation of the
corresponding patch on the temperature distribution surface in the

three-dimensional space where the third axis is temperature (see Figure

16).

A finite element grid diagram for the computer model mentioned in

Sec. 3.1 is shown in Figure 7.

In this diagram. the Y scale is exaggerated so that the elements in
each thin layer of the different materials can be shown including the

thin layer of glue. The reason to include the effect of the glue is its

qmmoz mmBszoo Mme mom QHmU Bzmzmqm meHth on mmDUHh

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

xi...
Amwvnhomv\\\\\\\
_ o S so 3 8 $ 3
38.28
38.35.7111 m . {Induce
433-381; J: 8 S m R m 12.?
33$..ng
unqmmmommfiwn on o: _: m:
Mmmflwmw .
I - 1 I mumm
83-2w: Ln: 0 N w ..
:13

/

mmmzsz A<Qoz mMmZDz Bzmzmgm

 

>—.¢_—————

35

low efficiency of conduction which might cause a significant temperature
difference between layers. The magnitude of thermal conductivity of
Devcon S-Minute epoxy is 0.389 IVE-K (data published by the manufactur-

er) compared with that of copper. 403 IVM-K. and 1.436 VIM-K for glass.

It can also be seen that the grid size becomes smaller near the
viewing hole on the grid diagram. because the accuracy of the tempera?
ture distribution within the viewing hole is more important than that

for the rest region of the system.

The comparison of the steady-state temperatures of T1. T,. T,. T‘.
and T, for the experimental results and the numerical results is shown
in Figure 8. All five of the operating conditions (see Table 1) are
included in Figure 8. and the parameter values used for predicting the
steady-state behavior of the computer model are listed in Tables 3 and

4. The relative error of Figure 8 is defined as

('r -'r )
Error($)= “'1 ”1 -1oos (3.1)

(Tmax. i-Tmin. i)

 

where T5,. is numerical result of Ti' and
Te.i is experimental result of Ti' and
T i’Imin.i is the total temperature range of the

III.

experimental result of Ti'

Figure 8 shows that all the relative errors (as defined above) are
within 105. and the relative errors even become smaller at the locations

of concern. T; and T3. So. the comparison of results reveals that the

36

 

RELATIVE ERROR CASE 1 CASE 2._.__.
1 CASE 3 ------ CASE 4—-—
15%
CASE 5 -—--—-'
10%

5%

 

0%

 

-10%
T1 T2 T3 Th , T5

45%

 

FIGURE 8. RESULT COMPARISON OF STEADY-STATE
TEMPERATURE DISTRIBUTION

37

steady-state behavior of the computer model with listed parameter values

can satisfactorily predict the behavior of the actual system.

3.2.2 PRCERAM FOR THE TRANSIENT TEMPERATURE DISTRIBUTION

This research places an emphasis upon gaining an improved under-
standing of the heat transfer characteristics of the cryomicroscope
conduction stage by making use of interactive computer graphics. Thus.
there was a need to present the various temperature distributions by
specific graphic outputs such as: three-dimensional temperature pro-
files. isotherm plots on two-dimensional graphs. and various displays of
the temperature histories at various locations within the cryomicroscope
conduction stage. No software was available to accomplish this task.

It was therefore developed by the author as part of this thesis.

The program deve10ped uses the forward-difference technique (Euler
method). because it is easy and simple for programming compared with the
finite element method for the transient case. This indeed greatly sim-
plifies the programming task and reduces the possibility of making
errors. One problem that does arise with this technique is numerical
instability. In order to avoid this problem. small time steps are often
required in the numerical integration of the governing equation. which

will significantly increase the computer time required.

As stated above a major objective of this research was to gain an
improved understanding of the heat transfer characteristics of a cryomi-

croscope conduction stage. The approach taken here was to examine the

38
thermal response characteristics of such a system using the interactive
computer graphic facilities (in this case the facilities of the Case
Center for Computer Aided Design of Michigan State University was used).
Since from the beginning of the present work it was clear that much
could be learned about the cryomicroscope system by considering a number
of possible design changes in the system. the software developed here

was designed to be very flexible.

It has been found that thirteen basic nodal equations can be
developed which will be sufficient to describe any rectilinear geometry
of interest ( i.e. any body constructed of orthogonal straight lines).
The program therefore has much more potential than has been tapped in
the present work. The thirteen equations also allow arbitrary boundary
conditions to be specified at any position. These equations are avail-
able as a 'library' available to the user to rapidly construct by

computer the heat transfer design of interest.

In summary. the program can handle these situations: any configu-
ration with rectilinear geometry. any magnitude of heat generation in
any region of this configuration. arbitrary boundary conditions. laminar

composite materials. and non-homogeneous grid size.

However. the primary purpose of this research was not aimed at
develOping software for a completely general heat transfer problem.
Consequently. the program still needs moderate modifications and proper
gridding for the numerical technique to handle all of these versatile

situations in an efficient manner. although some can be handled fairly

easily now. Extensive effort was not invested in improving the
pre-processing capabilities of this software at this time. A versatile
program is definitely necessary for this general kind of research
involving the effects of parameter variation and future effort is war-

ranted in this area.

The detailed derivation of the forward-difference equation for an
interior node will be presented in Appendix B. and the derivations of
the rest of the nodes will be omitted. Only the resulting equations are
presented in the following. and the corresponding nodal systems are

shown in Figures 9. 10. 11. and 12.

Forward-difference equation for an interior node:

(1+RKRY) [911.1, j-(1+1/ lune: fog“ ' jinx]

 

 

 

 

+
(1+nxmi';
RL’lef, 1-1-(1+RK/RY)OE' fare; 1+,Imn‘ +
H:
2gRL (1+RYRK/Ru)(91:31-01l )
= '3 'j (3.2)

Ailfilxguuxm 21m?

40

'-—dX1 +dX2 C—i

(I-I.J+1) (I,J+1) (I+1,J+1)

- dYZ ' Ei::>;>r- ' K2
I
Alf/A-

 

 

dY1 K1
(I-l,J-1) (I.J-1) (I+1,J-1)

FIGURE 9. INTERIOR NODE (WITH OR WITHOUT
AN INTERFACE AT THE NODE)

1+1

(A) NODB ON UPPER BOUNDARY

 

 

 

 

 

 

(I-l.J) (m) (lug)
K2
dXZ___J
-- dx-- --dx--
0 ‘ 7'4) 1
(I.J+1) K2 K2 (I,J+1)
;7' de 5;,-
(I+1sJ) (1-1 J)

(I.J) 4/5: Ba: ' V///r.(1.a)
Q \
EES§3 dY] E;:>L

K1 K1 '
(I‘D-1) l I (IJJ-I)

 

(B) NODE ON LEFT BOUNDARY

 

(c) NODE ON RIGHT BOUNDARY

 

dX1

1
(I

K)

dXZ

)
,J+1)

K2 dY

 

J\

4

(I-I:J)

$3

(1.

A C
J) (1*1.J)

(D) NODB ON mm BOUNDARY

FIGURE 10. BOUNDARY NODES

42

(I.J) (1+1.J) (I-I.J) (I.J)
¢(3 dY dY —
(I,J-1) ' (I,J-1)
dX dX
(A) NODE AT LEFT-TOP CORNER (B) NODE AT RIGHT-TOP CORNER

 

dX
(I,J+1) (I,J+1)
dY E23;7
(I.J) (I.J)
(1+1.J) (I-1,J)
(C) NODE AT LEFT-BOTTOM CORNER (D) NODE AT RIGHT-BOTTOM CORNER

FIGURE 11. CORNER NODES

43

(A) NODE AT LEFT-TOP INNER CORNER (B) NODE AT RIGHT-TOP INNER CORNER

 

(I,J+1) K2
de
(m) ,1 l
CI-I.J) QQEEEEE:;(I+1.J
& .L

K1

 

(C) NODE AT LEFT-BOTTOM INNER CORNER (D) NODE AT RIGHT-BOTTOM INNER CORNER

K1__dx1__f__dx2 Xm I dxa__i
(I. J+1) (I.J+1)
x2
dye dYa

\‘V
(I-1 J) 5‘ (I+1%

(I, J)

    

dY1

(La-1) K‘ ..L

EIGURE 12. INNER CORNER NODES

41+

Forward-difference equation for a node on the left boundary:

-(1+RKRY)(0'1'1-0n+11 1) -1111. (1+kYHe‘1‘11—02)

 

 

 

 

 

A31 mi.
RL'IO11-1-(1+ax/DY)O‘1'11+sref11+1lul 1 21m.
AY: K111151151
n+1
= (1+RYRx/RaHe111-o‘1‘11) 13 3)
211m?

Forward-difference equation for a node on the right boundary:

 

 

 

 

 

(1+RKRY)(01 1 1141111) -1111. (1+!!!) (9 111-0.. ) 1
Ail. ‘lfil
RL’le111-1-1(1+sx/BY)O11+axe111+1lnY1 1 2gRL
AI: x1'r1AY1Ai1
(1+RIRK/Rn) (91114-01111
= _ (3.4)
2mm

Forward-difference equation for a node on the upper boundary:

 

 

 

 

RL’(1+RKRX)(0111-10'1‘11) . m1n(1+u)(o1‘1 14)..) 1
4Y1 MAY.
[9111-11-(1+RK/RX)911+RK9111+1IRXJ 1 2312:.
Ali: x1T1AY1Ai1
(1+mx/Ra)(99f‘-e‘1‘ 11)
= (3.5)

 

2RtAt

1+5

Forward-difference equation for a node on the lower boundary:

 

 

 

 

 

:m.’u+mx> (01111-9111 1+1) 1. n1n(1+u)(91’11-e.) 1
AT: K1AY1
[011111-1-(1+nx/nx>e§11+no§1111/m 1 2sz
Ki: x1'r1AY1Ai1
n+1
= (“ml/Ru) (9111-01111) (3 .6)
21m?

Forward-difference equation for a node at the left-botttau corner:

 

 

 

 

 

 

-(91‘1 1-0'1‘1111) _ 31,3(0'1‘1 1-01‘1 1+1)
AX’ AY’
nL1(e‘1‘1 1-0...) - manna: 14)..) 1 2 gRL
m n? mum?
n+1
23m?

Forward-difference equation for a node at the left-tap corner:

-(91111 14911111) 11 21.49111 1_1-e‘1‘1 j)

 

 

 

 

 

Ai’ AI’
m‘x‘olil. 11.9“) _ HLxRL(01111 1-00) 11 2gRL
KAI m? snails?

«NW-(au )
= 1'1 i'j (3.8)
2RtAt

 

1+6

Forward-difference equation for a node at the right-bottom corner:

 

 

 

 

 

 

(91$1 11-91‘1 1) - 31,2(91‘1 1-0‘1‘1 1+1) __
Ai’ AY’
HL1(91‘1 1-9.) _ minus; 1-0.) 1 23m.
mi KAY H.113}?
«Wm-(3n )
= i" 1’ j (3 .9)
mm?

Forward-difference equation for a node at the right-top corner:

(911.11 1411.11 11) 1 1213(91‘1 1-1-0'1'11)

 

 

 

 

 

A'i’ AY’
111.1(9111 1-9.) - m.an(e'1‘1 re...) 1 23m.
mi n? u1AiAY
(91:343. 1’

 

= _ (3.10)
2RtAt

47

Forward-difference equation for a node at the left-bottom inner corner:

[2x11119111 11-(RKRY+mY/Bx+1/RX)91'1 1+(mY+1)e§+1 1 jI111]
A31

RL’ [Rxe‘1‘1 1-1-(RX+RKIRY+RXRKIRY)O'1'1 1+n(1+nx)o'1‘1 111/RY]

 

"'3
AY1

111.1(931 1-0.) _ m1RL(91‘11—e.) 1 23m.

K1Ai1 x1151 x1T,Ai1A§1

 

 

[RX+RYRK(1+RX)/Ru](0n+’-On )
= i’j 1” (3 11)
231A?

 

Forward-difference equation for a node at the left-top inner corner:

[911.1 11-(1+RKRY/ux+1/RX)91‘1 1+(RKRY+1)91‘+1 11m}

 

 

 

 

 

+
—3
All.
a n
m. I(1+n)e1’11-1-(RX+1+mx/RY)O'1111+mxe11 1+1/RYJ -
.17:
n n
K1A'i1 x1511 x1T,Ax1AY1

(1+RX+RYRKRXIR(1)(0'1""3-0111 j)
-- ’ ’ (3 .12)

 

22m?

 

1+8

Forward-difference equation for a node at the right-top inner corner:

l(1+RKRY)9‘1L111-(1+RKRY+l/RX)0'1’1 1+0?“ 11in]

 

 

 

 

 

 

A? +
1
RL'[(1+RX)0'1'1 1-1-(ax+1+nx/BY)O'1’1 1+sxe1‘1 111/RY] -
AY:
mxnnnefi 1-0.) shame: 1—0.) 23m.
- +
K1A'i1 x1A'i1 K1T1Ai1AT1
(1+RX+RYRK/Ra) (0'1”3-0’1‘ 1)
= ' ' (3 .13)

231A?
Forward-difference equation for a node at the right-bottom inner corner:

[(1+RKRY)01‘,11 1-(1+MY/n+m1)of1 1+mYe§+1 11121]

 

 

 

 

 

 

“1
sL’lo’1‘1 1-1-(RX/RY+1+RXRK/RI)0‘1’1 1+n<1+ax)e‘1‘11+1/RY1 _
AY:
manu(e'1'11-o,) _ 111.1(91'11-0.) 1 23m.
£11121 K1Ai'1 x1T,Ai1A'f1
I1+mr1u:(1+RX)llr.ul(9'1““3-(31I 1)
. ’ ' (3.14)

22111?

The forward difference grid for the computer model is shown in Fig-
ure 13. This computer program has several output features including:
numerical data indicating temperature of every node at each instant of

time; the rate of temperature change at T1 (see Figure 3 for location of

49

ammo: mmesmzoo mme mom QHmo movfimemHm Qismom

.m_ mmauHm

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

N. m m V m N _
x1411)
:9 mm m
.w _
m _
mnmmuy1. mums. “v _
zoEmm \r (“ABM “A w "
oneemmzmo m
e<mz mw
m
0_

 

50

T1): two- dimensional isotherm plots of the system; normalized
three-dimensional temperature distribution profile of the system: the
temperature histories of T1. T1. T1. T1. and T1. and the temperature

difference between T1 and T, (AT) as a function of time.

A typical output of numerical data for the temperature distribution
on a cross-sectional area of the cryomicroscope conduction stage at dif-
ferent time instants during the dynamic response is given in Figure 14.
In this figure. each data value is the nodal temperature at the cor—

responding position.

A typical output of the rate of temperature change of T1 also is
shown in Figure 14. The program calculates the average rate of tempera-
ture change during the interval between two subsequent temperature

distribution output.

Several representative two-dimensional isotherm plots at different
instants of time are shown in Figure 15. These graphs reveal how the
isotherms change at different time instants and the observer can get a
visual understanding of the temperature gradients and the direction of
heat flow within the system. This information can yield suggestions for
modifying the configuration as well as other parameters for designing a

better system.

Several plots of the three-dimensional temperature distribution
profile at different instants of time are overlapped and shown in Figure

16 and 17. On these graphs. the X? and Y-axes correspond to spatial

hmwdt‘dhfimhh
seedcmexmm
OOO-o—w-«anuun

\Tfifififidfidqdd

EDIIIIIIIIII
1")
m—mmnnmonou
msseeewemmm
«mmmmmwnnnn

"MMMMHHfic-‘d

éIIIIIIIIIII
LII
E0~Nflfih~hmm

f-I'dmI’II'I 1') J] 4253030:
1‘.) m :‘u‘ :‘d :‘u :‘d '1‘) In u‘: :n
L—mddvflflofivflddfi

QIIIIIIIIII

LI. OHID‘fo-IJ‘
H . . . . . .
"IPTI‘Q‘Q‘IJ‘
1‘; 1'! 1'! u‘w': In In
(1'13 daemon-eon
UM I I I I I I
1...
IT: «I‘d-(III‘IOO
('3 . . . . . .
~— I'r-OIJZIOM-e
a: I'II'III‘MII'Oeu
'- "Me-1H0!“
'J'I I I I I I I
H
Q mm
m iv
a“; «dell
1:. 040-4
In. I I
'I
"1 9303
Lu . .
I1, :I‘HI‘)
r- '0‘!)
m use-s
— II
N
.L

IS

Hume

Q
J

TEJfERATUHE

RAIE OF

51

mmmmma
M~§I~§I<I~lr§
wmmmmm
I I I I I I
mmmmmm
ITICIIIVIITIITIIN
mmmmmm
I I I I I I
nflflflfld
mammmm
000900
IIIIII
¢¢nnnm
.(‘Idx‘jd-I
)OOOOCJ
.fl

I

’I

.—1.',

Odfd‘TDOf'JO‘CDEU
deedcmesmm
«undue—MIDIDIDHT
flunndddnmdd

«IIIIIIIIII
9
nmemnmo~¢oo
mwswseesmmm
«mmwmmmmmmm
Mdmdfiddfidfl
"IIIIIIIIII
LI
Emmemoommee
hflflflflfiflflmfifi
isfcrcrufsfuuiuvufiui
t—v-eo-so-so-nouo-qo-qo-eo-ao-e

’«IIIIIIIIIIII

2;) N n: ~t In 0 O
H . . . . .
"I: "f (‘0‘th
L“ 1') fun :I'. 43.13
D «Manda-soc
M 1 I I I I I
1..
D I‘df'lh'fldsm
M IT-IZIU’JCJ-u-u
m mnnnno
3" “H"No-H
III I I I I I I
I4
G 00
m 14
g flfl
-1 «as
r- I :
'I‘
a; :!'-IJ~
LU . .
.‘11 WM“
p... Jpn
“J «so-s
~ II
I..1I
.L

.1 C/SEC

“J
fi

"0.

IS

H-’-.NCE

I

g
‘

123,315 0? TE:'!L-‘ERATIJ;~£E

nmmmon
dddddé
000000

undo-1'1"“

I I I I I I

mennmmnooo
cacddnemmm
o-I—I—I—I—Iv-HI‘IDDII‘I

mddflddddddd

~‘I‘I I I I I I I I I I
‘3'

CIUIQP-D‘Odl‘dD—I-I
nesssmesmmm
on I‘d I‘d W I“) I1! 1") In In In In

"ddfiw‘flHHO—‘d
:IIIIIIIIIII
II
.1
ZL‘JQU‘IDI‘JOFIOHJ‘ID‘I
k-I'II‘II'II'Irrr-sflmthf"

“."J :‘d If: :‘I.’ :‘d ('1: u‘. '1‘) In In
:mmmflddcqo-soqdd

'IIIIIIIIIIII
VIVID-Odd

9‘.) ~14 r; I)? (j (73
2'31‘32‘1 In as] «J
as -s .4 as as as

Iii:

memes;
NHNflflfl
HMO-Q0404"

IIIIII
00
II‘IIIT
0'0

on!”

I I
QC)
J) «‘i
.13 .13
'1."

TErff'ERATUFiE D I S I P. 1 DUI ICéI‘I

IHE

Inwmmmru‘

I?
W

'3

1‘3

CHf-IIIJ E

RATE OF TEnFERATUHE

FIGURE 1h. NUMERICAL RESULT OF NODAL TEMPERATURE DISTRIBUTION

C mm<ov 95qu :mmmaomH A<onmzmzHDIOBB om— mmeHh

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Q 3m 20:: H0 9% a 0::
4m
0 r I w -
03m 02m;
.3 .mm .8
2
R/
Hm 9.2 20:: “mm 3 m. 0.1:
.a
v _ _
.mmé .8 .a

Am Mmdov ZOHBDmHmBmHQ mmaedmmmzme d<onmzmszlmmmmB QMNHQ¢EmOz owp mmDUHb

o." m.-

uww m.~_ m~ mx~p m.-

me¢qa mammoo

000: ozHamH> m .

 

7)
r)

 

 

 

 

 

 

 

 

 

 

 

 

. Aaz<m xm<mv
mme<mm mm<0o.:////1: onomm .. ...
one¢m020c e<mm

 

ho ZOHB<DOA

 

mqmz<m 0<0Hooqum mo one<ooq ;

1m»

Am mm<ov mez<smzH mzHe ezmmmmmHm e< monesmHmemHm
ampe¢mmmzme 0¢onmzmsznmmmme nmNHg<zmoz mo oneHmommmmsm .m— 3:230;

9.. m.s

II‘ b F

x I. I . b

m s

 

umm a.m© m. m:.»

 

 

 

 

 

 

.ozn“ so "0

 

 

 

 

 

5..

 

 

 

xmh

55

axes in Figure 13 and the third axis corresponds to a normalized tem-
perature scale. Since all three axes are normalized. the observer can
get a feeling of the relative magnitude of the effect which a given
design change may have on the temperature distribution within the sys-

tem.

An exciting point to be made is that these graphs of different time
instants can be displayed in a dynamic fashion on certain high speed
terminals such as the Tektronix 4114. This equipment can be used to
produce a computer-generated movie which presents a vivid and dynamic
impression of the transient response of the temperature distribution by
using the DRAM.A software. which is available in the Case Center of
Michigan State University. Preliminary work of this type was accom-

plished as part of this thesis.

A typical output of the temperature history at a specific location
is shown in Figure 18. The major purpose of this type of output is to
compare with the experimental results of temperature histories at the
corresponding position in the actual system (see Figure 5). The rate of
temperature can also be calculated from the slape of the curve on the

chart. if necessary.

A typical chart of the temperature difference. AT, as a function of
time is shown in Figure 19. This chart gives the user a convenient mean
of inspecting the influence of parameter variations on the important
quantity. AT. which is an indication of the uniformity of the tempera-

ture distribution within the viewing hole.

56

A_ mm<ov _e mo emoemH: mmae<mmmzme .w_ mmauHm

..~m_ -.mm_ p.00. o.-_ o.wm o.~n -.m¢ ..*N c.-

.uwm.

 

T

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

mxuh m>

d

mm3p<mmaxmp

o..m
u.m~
..-m
..mm
c...
..m<
..-m
..mm
-.-o
a.mo
u..~

.u mmawmo.

57

A. mmeuv mzHe mo oneozsm
e me me oz< _e zmmaemm 0020mmmmHm mmseemmmzme .m_ mmouHm.

.uwm.
..~m_ ..ow_ -.¢<_ ...~_ 9.0m ..~n ..o¢ ..¢~ o.-

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

mxak m> muzmmuuu_o xmh

3..
.5.
.5.
3.
..~.
....
..~
3
.6
.5
..:

.u mmmwmo.

CHAPTER 4

COMPARISON AND DISCUSSION OF THE EXPERIMENTAL AND NUMERICAL RESULTS

4.1 COMPARISON OF TEEyEXPERIMENTAL AND NUMERICAL RESULTS

Comparisons of temperature histories at different locations within
the cryomicroscope conduction stage for each case (see Table 1) are
shown in Figures 20-44. In these comparisons. the experimental tempera-
ture histories are presented with dashed lines. and the
computer-predicted temperature histories are presented with solid lines.
By direct superpositioning of these two curves, the deviation of the
numerical results from the experimental results can be visualized

immediately.

A quantitative evaluation of this deviation is presented in Table

5. The maximun relative error (M.R.E. in Table 5) is defined as

58

59

 

"943mmm HdonMZDZ

”Egsmmm Adezmszmmxm

. mm¢o mo Apv.e .om mmmon

..mm_ ..ow_ ..<<_

.umw.

..-_ ..mm

..~n -.m< -.em u.-

 

 

 

 

 

 

 

\\
‘Ux

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

mxah m>

 

 

 

 

 

 

 

ma3h<mmaxwp

..o~
..m~
-.-n
..mm
...e
..m<
..-m
..mn
..-o
o.mo
...n

.u wwmwwo.

60

 

”Bdsmmm A¢DHmmzsz

“Egbmmm adezmszmmxm

_ am<o mo Aevme ._m mmaon

.uww.
—.mm_ —.oo_ -.<<u ...~u -.mm -.~s ..m< ..<~ u.-

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

wxuh m> m wank<mwanF

o..~
..m~
..-n
-.mm
...¢
-.m<
...m
..mm
...o
..mm
...u

.u wuamwo.

61

. mm<o mo “pvme .mm mmauHm

8mm.
..~m_ .60. .4: .63 ..8 -.-. 90¢ .JN o.-

 

gin
\ ..m~
\\\ II]. ..-n
.1... :1...qu ..mm
2:
29 8 H68.
“eqammm m<onmzaz ...m
, .6»
“530m 05.20szQO 28

 

 

 

 

 

 

 

 

 

 

 

 

 

.60
2:

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

ms: m> a $5?ng

62

_ mmeo mo “00:9 .nm mmauHm

8mm.
..~m_ .62 .4: ...~— 98 ..~s 99 .4m .6
_ _ _ I -
III... I ll ..mm
98
..mn
9:

 

 

 

 

 

 

 

 

fine 8 $8.
“902mmm 0<onmzsz ..-m

 

 

. .3
H mm“. 4 J 9 NW «N 0
ea: c m¢ezg20mm x: . .m

 

 

3mm
22.

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

ms: m> 0 9.5.25ng

_ mm¢o ho Aevme .:m mmsuHm

.uwm.
o.mo_ ..mm_ ..¢<_ ...~d ..om ..~h ..m< ..<~ o.-

 

...~
3m
3m
|.. 3”
k.\\ ..:
I 3. a was.
”Baammm 64055222 ..-m

 

 

 

 

 

 

 

 

 

 

 

 

lllll .u ..mm
”53mm 05.20280- ”a 2.0

 

 

9mm
.65

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

U3... m> m 952185»

6h

 

”BADnmm Adonmzaz

“Eqsmmm adezmszmmxm

m mm¢o mo ApVPe .mm mmson

.uwm.

..~m_ ..mm_ ..<<_ ..-d ..mm

..~n ..o<

..¢~ o.-

 

 

 

 

 

 

‘3‘
‘b

 

 

 

 

 

‘R-

“‘A

.\

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

wrap m>

d

umbp<mwaxmh

o..~
..mm
...m
-.nm
...q
..n<
..-m
-.mm
o..o
-.mo
...s

.u muzwwo.

65

 

"Eqsmmm A<0Hmm232

“Egbmmm A¢Bzmszmmxm

m mm<o mo Apvme .mm mmauHm

-.~m“ -.mm_ -.<v— -.-¢ -.mm

.umw.

o.~n ..o< ..<~ o.-

 

_

 

_\

 

 

‘I

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

mx_p m>

m ma3p<mmaxwp

-.-
..mm
...m
o.mm
o..<
..me
...m
..mm
..-0
-.mm
..-s

.u mwmwmo.

66

 

"BADmHE A¢0HmMZDz

"BADwmm Q<Bzm2HmmLxm

m mm¢o mo Apvme .um mmauHm

-.~m_ ..mmu -.¢<. ...~_

.umw.
-.mm

-.Nh ..o¢ ..e~ ...

 

Q

 

 

 

Ulllv‘u

 

 

 

 

‘

\

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

mxHh m>

m mmap<mmaxmp

..-
..mm
...m
..mm
o..<
..m¢
a..m
o.mm
...o
..mo
...s

.u mmawwo.

67

 

"BADmmm A<0HmMZDz

”Hg

0mm A¢BzmszmLXE

m mw«o mo Anvse .wm mmpuHm

..~m¢ ..an_ ..<<_ ..-_ ..om o.~n ..m< -.<~ c.-

.umm.

 

A

 

 

 

III
I ,

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

m:_h m>

d wmap<mwaxwp

..-
..mm
...m
..mm
..-<
..m<
...m
..mm
...o
..mm
...s

.u wuawwo.

68

 

"945mmm A40Hmmzaz

"Edammm A<Bzm2Hmmuxm

m mm<o mo Anvme .mm amouHm

.uwm.

-.~m_ ..mm_ o.<<_ ...~. o.wm

-.~h

-.o¢ ..<~ c.-

 

# .

 

 

 

 

j§
I

 

 

 

 

 

 

 

 

 

 

\1

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

mxwh m>

m wanb<mmaxwp

o..~
..m~
..-m
..mm
...<
-.m<
...m
..mm
-.-0
..mm
...u

.u mumwwo.

69

 

"Edammm Adonmzzz

"Eqsmmm A¢Bzmszmmxm

m mm<o mo AeVFB .om mmauHm

.uww.

-.~md ..mo_ ..<<_ ..o~d ..mm ..~w ..m< ..<~ ...

 

 

 

 

I'IV'J "]

 

 

 

/I

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

mxuh m>

.

maap<mwasmb

....~.
..mn".
..-m~.
..m~_.
...-..
..ms. .0 mumwwo.
...m.
..mm.
o.-
..mm
-.-m

70

m mm¢o mo Apvme .Fm mmauHm

.uww.
..~m_ -.om_ -.<<— ..-_ -.om o.~n o.m¢ -.*~ g.-

 

..-~.
-.ms_.
.62.
..m~_.
......
..ms. .0 mmzwwo.
"SHIRE $on0222 93..
IIIIIII [1, ..mm.
”aqsmmm 0<ezmszmuxm ...
..mm
...m

 

 

 

 

III '

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

ms.» w> m maah<mwaxmh

n mm¢u mo Asvme .mn mmame

.89
Quad .69 .43 .63 98 9mm .dv .Jm c.-

 

.63.
.62.
I 1.2-
.63.
UK 9.3.
x 9mm. 8 EH:
”emsmmm 0<onmz=z r, ...m.

’7 ..mm.

"903mmm 0¢Bzmszmmxm ...

 

 

 

 

 

,9

 

71

 

 

JP

 

 

 

 

 

 

 

 

 

 

 

 

 

9mm
_ 92».

 

 

 

 

 

 

 

 

 

Us... m> m 595%....me

72

m nm¢o mo ApVJB .mm mmauHm

.uwm.
..~m_ ..auu ..<¢_ -.-_ ..mm o.~s ..m< ..e~ ...

 

....~-
..ms_.
..om~.
/ / 9mm?
x / .62.

 

 

 

 

 

 

 

 

I

I ..mu. .0 Hyman.
"SEE 053,052 x, 93..

..mm.

"905mmm q<azmszmuxm ...

..mm

...m

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

wx_h m> e wm3h<¢maxmh

73

 

”eqommm A<0Hmmzoz

"Habmmm A¢Bzm2Hmmaxm

mm¢o mo Apvme .am mmpon

.uwm.

-.~m— ..owd ..<<_ o..~. -.om ..~h

..m< ..g~ o.-

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

I”
IIIIIII.:1.J1I I]/{
III)
[‘1’ /
,0
r I
I,
wink m> m ma:»<awaxmp

..-~.
-.ns—.
...md.
..mmd.
...-E.
..ms.
...m.
..mm.
o.-
..mm
...m

.u mummuo.

74

: ammo mo Aev_e .mm mmpuHm

.uwm.
..~md ..mwu ..<¢u ..-d -.mm -.~h ..o< -.¢~ o.-

 

..-.~.
..msu.
...mfi.
..mm..
is} .52.
/ 92. a .058.
“003mmm 0<onmzaz I ...m.
IIIIII I: .1 ..mm.
"902m0m 0<ezmszmuxm , ...
..m~
2....

 

 

 

 

 

 

 

 

 

I

 

11’

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

ms.» m> d mm3p<amawa

75

: xm<o mo Auvme .mm mmame

.uwm.
..~m~ -.aa_ -.<<¢ -.-_ -.mm ..~s -.o¢ -.¢~ ...

 

-.-~.
..mnd.
..gmd.
..nm_.

:flffl ...:.

 

 

 

 

 

II 1.1.11 My

 

 

 

 

. —.m5. .0 muuwwo.
”emsmmm 040Hmmzsz ...m.

I IIIII I I. 3....

"eqswmm 0<ezmszmnxm I ...

— ..m~

- ...m

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

mx_» w> m wm3h<mmaxwp

76

J 00.0 mo ...ns .um mmson

.HH:
..~m. ..oo. ..<¢. ...~. -.mm ..- ..m< -.¢~ -..

 

-.—Q~.

..ms.-

...m..

11114111111] .62.

J74” ...:-

711 ..ms. .uumwmmo.
“eqammm 0<0Hmmzaz ...n.

lllll I I_ ..m~.

"emammm 0<ezmszmuxm . ...

..mm

23

 

 

 

 

 

I

 

 

 

I

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

mxah w> m w¢35<mmnzmh

77

: mm<o mo .ev:e .mm mmson

.uww_
-.~m_ -.m®_ ..e<. ...m. ..mm o.mn -.m< -.<~ o.-

 

_ ....~.
-.m~..
I .62.
I .53.
, ,/ ...:.
Aw TI; ..m~. .0 Human.

"52mm... 0335.52 .2 ...m.
m/ .;m~.
“egpmmm 0<eszHmmmxm , ...
..mm
...m

 

 

 

 

 

I},
I

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

us.» w> a wm3h<mmaxwp

78

 

"BA mum A¢0Hmnzaz

"8Q mmm Q¢Bzmszmuxm

-.~m. ..mo— ..¢<— ..o~_ ..mm

: mm<o mo Auvme .mm mmpon

.uww.

..~h ..m< o.¢~ ...

 

 

 

 

 

 

U‘

 

I

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

mxub m>

m wmob<amaxwp

..-o~-
-.mn—.
...m.-
..mm..
...-..
..mn. .u mmawwa.
...m.
..mm.
c.-
..m~
...n

m mm<u mo A...e .o: mmpmHm

8mm.
.60— ..8— 6:3. .62 6.8 .65 .60 .JN .6

 

.63.
.62.
.62.
.63.

 

 

 

 

79

illwnpl.... ....—.
nnurlj I, ..mu. .0 mmawmo.
"Sammm Roamzsz I .. .m.

IIIIII /f ..m~e
"egammm 0<Bzmszmnxm / ...

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

.6m
_ 1

 

ms: w> . £523.!“—

80

 

”Bu mmm A<0Hmmzaz

"BA

mlq‘d
IIF

m Adazmszmmxm

m mm¢o mo Agvme ._: mmsuHm

.uwm.

..mm. ..oo. ....d ...m. ..om ..ms ..o. ..e~ ...

 

 

 

 

 

 

 

Alli

 

I
I
%/1
I

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

mxah m>

m wm3p<mwnxmp

....m-
..m~..
...md.
..mmd-
......
..ms. .0 mumwwo.
...m.
..mm.
...
..mm
...m

81

 

”9420mm A¢0Hmmz=z

’. .’ I I ... '

“Basmmm Q¢Bzmszmuxu

m mm<o mo Auvme .m: mmsuHm

..~m_ ..mw. ..vv— ...~_ ..mm ..wn ..o< ...m ...

.uwm.

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

ms.» w>

m wm3h<mmaxwh

....N-
..m~_.
...m_.
..m~_.
......
..mu. .u wwmwwo.
...m.
..mm.
...
..mm
...m

82

m mm<o mo Apvze .m: mmauHm

Bum.
..N2 .69 ...: ...m. ..8 ..mn ..mv ..dw ...

 

_ _ _
.l...l ......a _ ..m:.
L ....-. iflyr. ...2.

40/ .....2.
x 2 ...3.
/ ..fi. a E.
...m-
..a.
...
..a
...m

 

 

 

 

 

 

 

 

I
d

 

“Edammm A¢0Hmmz=z

 

 

"BADMQm Adezmszmmxm

,49/1

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Nah m> d map—(mung.

83

m mm<u no “gyms .s: mmson

.umm.
..~m_ ..mo— ..e¢_ ...m. ..mm ..Nn ..m< ..em ...

 

...-N.

..m~_.

.....-

..m~_.

......

I 1.74::/ ..mn. 8 $58.
"53mm 203352 /7 ...m.

.. :- :....u x. ..mm.

...

_ 1,114“; ..wm

 

 

 

 

 

 

 

 

 

 

 

"

"Eqawmm A¢Bzmszmmxm

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

...m
mx_p m> m wmah<mwawa

84

TABLE 5

QUANTITATIVE EVALUATION OF RESULTS COMPARISON

Nggimnm Rglatize Eggor §l,R,E,2

T1 T: T, T‘ T, Average
Case 1 4.7% 3.8% 7.4% 30% 3.1% 9.8%
Case 2 5.0% 4.5% 7.1% 23.7% 5.6% 9.2%
Case 3 5.6% 5.1% 9.9% 14.9% 14.0% 9.9%
Case 4 8.6% 4.3% 11.0% 11.7% 12.5% 9.6%
Case 5 12.9% 8.6% 10.3% 16.3% 11.4% 11.9%
Average 7.4% 5.3% 9.1% 19.3% 9.3% 10.1%

85

(T -T )
n,i ,1
M.R.E.(%) = ° "x ~1oos (4.1)

Imax,i'Tmin,i

 

where Th'i is the computer-predicted temperature history of Ti'
and Te,i is the experimental temperature history of Ti'
and the subscript ”max" in the numerator indicates the maximun
difference between the experimental temperature history, T;'i(t).
and the computer-predicted temperature history. Th.i(t).
and Tmax,i is the highest temperature reached by the
experimental result of Ti'
and Tmin.i is the lowest temperature reached by the

experimental result of Ti‘

4.2 DISCUSSION

4.2.1 DISCUSSION OF RESULTS COMPARISON

Generally, most of the maximun relative errors are within 10%. The
major concern of this research is whether the computer model can predict
the "trend" of the behavior of the actual system rather than specifying
an absolutely accurate prediction for the temperature distribution withr
in the cryomicroscOpe conduction stage. A more accurate prediction is
possible. but it will cost extra expense due to introducing precision
measurement facilities. detailed modification of the computer
model....etc. So. according to the present purpose of this research

this amount of deviation stated above is considered as acceptable.

86

In cases 1 and 2. the system has been heated only (see Table 1).
and the numerical results are quite satisfactory. because most of the
maximum relative errors are below 7.4% except T‘ (.g. Fi‘urg 3 for the
location Of T4). Apparently. the results fail to match at T} where the
maximun relative errors are 30% for case 1 and 23.7% for case 2.
However. this amount of deviation represents only approximately 3°C
because of the small total temperature range of T;. Fortunately, T}.
the temperature of the capper plate. is not of major concern. because
the biological samples are placed on top of the cover glasses. It is
the temperature distribution within the viewing hole that is important

and the errors in this vicinity are below 5%.

In case 3. the system has been cooled only. and the experimental
and numerical temperature histories match quite well. within 7.5%. A
comparison 0f thO 1°831t8 it T4 still shows that the model is inaccurate
here. but it is improved in this case (14.9%) as compare with the 30%

and 23.7% observed in cases 1 and 2.

In cases 4 and 5. the system has been cooled and heated simultane-
ously. These cases were designed to serve as a check on the predictions
of the computer model as constructed using cases 1. 2. and 3 as compari-

sons where only one kind of stimulation has been imposed on the system.

From Table 5. the error for T1 of case 5 is 12.5%. This error
porbably arises from a singular error in the experimental data.
Returning to the chart of results comparison (Figure 40). and it can be

found that the experimental curve for this case is not as smooth as for

87

other cases where the errors for T1 .1, much gngllor.

This situation is relevant to one of the possible weakness in the
present research. Only one set of experimental data was obtained for
each case. Multiple experiments would be an improvement by increasing
the statistical basis for describing an average experimental response in

each case.

4.2.2 ERROR DISCUSSION

This section summarizes some potential reasons for why the relative
errors may occur even after parameter variation is allowed in the com-

puter program in order to minimize such error.

One reason that such errors may exist is that a large temperature
range (200°C) is distributed over a small system (a few centimeters).
making temperature measurement errors more likely unless precision meas-

urement techniques was employed during the experiment.

Manual decoding of the analog signal of the thermocouple from the
strip chart into the temperature scale may also cause some experimental
errors. Since the finest division of the strip chart approximately
equals to 2.5'C when it is converted into temperature scale with an
Operate temperature range from -200°C to 50°C, the inaccuracy of the

manual decoding is t 1.25°C.

Another possible experimental error arised during the conversion of

88

analog signal (millivolt) into the temperature scale ('C) is related to
the thermocouples and the conversion table (reference [11]) used.
Generally. the thermocouples (Copper-Constantan) used in the experiment
should be calibrated according to the conversion table. so that the
conversion table can be applied to the conversion of the analog signal
of the thermocouple into the temperature scale. No calibration of the
thermocouples used in the experiment has been done. so this becomes

another source of experimental error.

Errors may also arise from the perturbation due to the act of mak-
ing the measurement. For example; to measure the temperature of the
viewing hole. a thermocouple has to be put on top of the viewing hole.
Since the width of the viewing hole is only 0.76 cm. there are several
concerns. One concern is the accurate placement of the thermocouple.
Another concern is that the epoxy used to fix the thermocouples in place
may change the local heat transfer. In addition. the thermocouple
itself acts as a fin. Thus. a difference between the actual temperature
of the position.where the temperature is intended to be measured and the
temperature at the junction of the thermocouple may arise leading to

errors.

When the system was being cooled. it is assumed in the computer
model that a step input of the cooling occurs (that is an instantaneous
change in ambient temperature and corresponding heat transfer coeffi-
cient. H). but this is not the case. During the experiment. the valve
was Opened and coolant was introduced into the system. The coolant has

to cool the channel down to close to the coolant temperature first. It

89

also takes time for the coolant to go through the whole system due to
the finite flowrate. After these transient events the coolant can then
perform the expected cooling effect on the system. The actual case will
show some delay in the cooling effect. and this can be accounted for by
a slightly higher temperature at the beginning on certain experimental

curves in comparison with the numerical curves.

One of basic assumptions mentioned before is that the temperature
gradients along the Z-axis are negligible (see Figure 2). However this
system is only 2.5 cm in this direction. and this will certainly influ-
ence the result to some extent. Furthermore. the present situation
where a large temperature range is distributed over a small system tends

to magnify this deviation.

Due to the turbulent nature of the coolant flow. which is a
two-phase flow in this case. the flowrate is quite unstable. This in
turn will introduce some fluctuation into the experimental data because
of the fluctuating of cooling effect caused by the fluctuating coolant
flow. This fluctuation will be damped slightly by the intermidiate
material between the coolant channel and the locations where the system

temperatures are measured.

The errors due to the forward-difference numerical technique may
also contribute to the deviation between the experimental curves and
numerical curves. The magnitude of this error depends on the time step
taken and the distribution of the nodes selected for the

forward-difference equations.

90

It is also assumed that the thermal prOperties during the operating
temperature range remain constant. but they actually vary with tempera-

ture, and this assumption surely has some contribution to the errors.

Heat generation in the actual system occurs within a thin film of
metallic oxide coating. which can be treated as a line on the
cross-sectional area diagram (see Figure 3). But in the computer model.
the heat generation region is treated as distributed over a band com-
posed by a series of rectangular regions each centered upon a node (see
Figure 13). This results from the manner in which the heat generation
is introduced into the nodal equations for the forward-difference method

(see Appendix B).

Furthermore. some thermal pr0perties of the material used to build
the system are not available. The values gained from trial and error on
the interactive computer terminal are not likely to be the actual values
for these thermal prOperties. This situation also will cause some

errors.
After reviewing these potential sources of errors. it is considered

satisfactory with the computer model that the numerical data is within

10% deviation of the experimental data.

4.2 .3 DISCUSSION OF THE GRAEICAL (IJTPUTS

The graphical outputs of the normalized three-dimensional tempera-

ture distribution profiles and the plots of isotherms yield valuable

91

information. For example:

In cases 1 and 2, the system has been heated only. The isotherm
plots and the normalized three-dimensional temperature distribution pro-
files at several instants are presented in Figures 45 and 46. In
general. the temperature of the entire system is increased by the heat
generation in the glass heater. The temperature distributions of glass
heater and cover glasses are apparently tilted» and this suggests that a

significant temperature gradient exists in this region.

On the contrary. the capper plate has a fairly uniform temperature.
although the temperature is slightly increased in the right side of the
plate closer to the heat generation region. This temperature uniformity

is due to the high conductivity of the cOpper plate.

There is a significant temprature drop (about 10°C) across the gap
between the glass heater and the Kapton film (see Figure 3). This
emphasizes the necessity of taking care of the contact resistance

between the glass heater and the Kapton film.

The temperatures around the viewing hole are higher than the tem-
peratures of other portions of this system. and this indicates that the
thermal resistance of the air is much greater than that of the contact
resistance with the Kapton film. because the lower surface of the view-
ing hole is exposed to air while other portions are in contact with the
Kapton film (see Figure 3). This high thermal resistance of air leads

to the accumulation of internal energy around the viewing hole region.

92

Am 9.: _ mmmfi: mammEomH ..Ho mmzommmm 5&3?ng .m: $.5on

gum m;m mawx:

 

 

 

 

 

 

 

. Fe

.mm.m~
.3 .mm

 

umm ..9 2me

4m

 

 

 

 

 

 

 

 

8:

.mm .«m
.mm

 

umm ...‘m 8me

4m

 

 

 

 

 

_.

.mmdmém

 

 

 

 

 

 

.mm

umm ..m 8 ms:

 

 

 

 

 

 

 

 

93

AN 024. _ mmmddv ZOHBDmHmbHQ mmDB/Nmmmgme le ho mmzommmm BzmHmH/Eme so: mmeHh

.._ m..

 

mum ..mw m. m:_p

 

 

 

 

 

 

 

 

 

:w»

91+

and the temperature around here becomes higher than in other regions.

In case 3, the system is cooled by the coolant only. and the
three-dimensional graphs and isotherm plots are shown in Figures 47 and
48. Basically. the temperature of the system is reduced by the cooling
effect of the coolant such that thermal energy diffuses to the boundary
contacting the coolant channel with the rest of the system modelled.
When the cooling effect is penetrating to the center region of the upper
half of the system (glass heater and cover glasses). the temperature
distribution within this upper region becomes V-shaped. because natural
convection on the left boundary (see Figure 6 and lower surface of the
viewing hole (which is on the right portion of the computer model) act

as heat sources in this case.

The temperature around the viewing hole is still the highest tem-
perature region of the system and this means that the heat flux comes
primarily from the viewing hole region and flows toward the portion con-
tacting with the Kapton film. Passing through the Kapton film. the heat
flux continues to flow down to the copper plate. and travel to the left

where the coolant channel is located eventually.

This description also can be inferred from the isotherm plots which
explicitly display the temperature gradients. The temperature gradients
reveal the heat flux within the system since the flux lines will be nor-

mal to the isotherms.

In cases 4 and 5. the sytem is heated and cooled simultaneously.

Am mm<ov mammzeomH ho

umm ..mm 2 mi:

.53.

mmzommmm BzmHmz<mB

 

 

 

 

 

 

 

 

 

Ker

.mmém-

r7
0..

umm 9.9 m. m::

.mn-

 

 

 

 

 

 

 

/ vi

.mm-

on: mmeHm

uwm ..mm 3 m1:

.mm- ..3.

/L1

.3.
.mm.

 

 

 

 

 

 

 

 

umw ..m 8me

.mm. .3.

 

 

 

 

 

 

 

 

96

Am mmeov oneDmHmemHm ampeemmmzme mum ho mmzommmm ezmHmz<me .w: mmsuHm

... m.

>
m: C \ \m.

... m.-

 

 

 

 

 

 

 

 

.omm :m "all

 

 

 

 

 

.85. N up

5.

97

and the normalized three-dimensional graphs and the isotherm plots are
shown in Figures 49 and 50. The basic configuration of the temperature
distribution is similar to that of case 3. except that the high tempera-
tures around the viewing hole are magnified further due to the

superposition of the heat generation effect on the previous results.

At the beginning, the temperatures within the viewing hole region
will be somewhat greater than the ambient temperature. However. eventu-
ally the temperatures will be reduced by the cooling effect. Also. the
temperature gradient becomes steeper and this is displayed on the isoth-
erm plots by a higher density of local isotherms. and a steeper
temperature distribution profile on the upper portion of the computer

model on the normalized three-dimensional graphs.

The normalized three-dimensional temperature distribution graphs of
all the five cases have one similar aspect. That is. the temperature
distribution within the viewing hole is a parabolic curve. This can be
illustrated by Figure 51 which refers the steady-state nodal tempera-
tures within the viewing hole of all the five cases. and plots them on
one figure with pr0per normalization as indicated in Figure 51.
Although only the steady-state nodal temperatures are plotted. the
dynamic nodal temperatures within the viewing hole region display the

same parabolic characteristics.

A curve of the parabolic equation

Y = 1 - X’ (4.2)

98

Am 92¢ : mmm¢0v mzmmmaomH ho mmzommmm BzmHmz¢mB

umm ..mm 8me

.mm..

 

 

 

 

 

 

 

 

:1

. . .ms.
...nm-..m.

 

 

uwm ..m: mums:

.ms.

 

 

 

 

 

 

 

 

om: mmstm

umm ..«m flux:

...?

 

 

 

 

 

 

 

 

 

1W

:

umm ...: 8me

..m.

 

 

 

 

 

 

 

LLII

.8 ..

 

Am 92¢ : mmm<ov onesmHmemHa mmpe<mmszme mum mo mmzommmm ezmHmz<ma .om amouHa

 

 

 

 

 

 

 

 

 

... m..
.
cum ..mw m. ms.» 9
E
... m..
I‘I .

’\“

. 1% I“
...

somw N H#\

twp

mqom wzHBMH> mma szBHs ZOHBDmHmBmHQ mmae¢mmmzme UHAom¢m<m Sm mmwah

mqom wzHBmH> mo mmasz mme 20mm w¢3¢ muz<BmHQ QmNHA<fimoz

100

 

 

 

 

.O N.O 3.0 0.0 moo .—
mqom UZHBWH> Nvfl
N.O
\
J \
00 a
2 2 \
m... T \
. v. 0.0 \.
mu. m...
__ \\
Y
.B\
\ \ mx .. P .1 mmfi
\B\ n mmg
\\. N mm<o

OU<1>O

101

is also plotted on Figure 51 for comparison. It appears that the tem-
perature distributions within the viewing hole of all the five cases
almost exactly match the parabolic curve. This also matches the results
of reference [4] which has some discussions about the parabolic charac-

teristics of the temperature distribution within the viewing hole.

In summary. although the abolute values of temperatures in the
present model may be somewhat in error. this is considered acceptable.
Having established a computer model accurate to within approximately 10%
or better when compared to experimentally measured data. the major use
of the model will be to predict general trends to be expected of the

thermal response of the system due to possible design changes.

CHAPTER 5
EFFECTS OF PARALETER VARIATIONS ON THE HERMAL RESPONSE OF

A CRYOMICROSOOPE CONDUCTION STAGE

5.1 RESULTS OF PARAMETER VARIATIONS
There are seven possible design change areas to be considered in this
section of parameter variations. They are:
1. Width of Viewing Hole 2. Thickness of Glass Heater
3. Magnitude of Heat Generation 4. Width of Heat Generation Region
5. Thickness of Copper Plate 6. Material of Metal Plate

7. Position of Coolant Channel

Each design change will be evaluated according to its effect on:
(a) the temperature difference between T1 and Ta (AT), which indicates
the uniformity of the temperature distribution within the viewing hole

region; (b) the maximum heating rate at T1 where the biological samples

1023

103
are placed and; (c) the maximum cooling rate at T1. These three quanti-
ties. AT. maximum heating rate. and maximum cooling rate. are the design

criteria for the cryomicroscope conduction stage (see Figure 52).

During the investigation of parameter variations. most of the
parameters are still the same as those which appear in Figure 6 and
Tables 3 and 4 for case 5. but a particular heating and cooling process
is applied on the system. The initial temperature distribution is set
to be the same as the steady- state temperature distribution of case 3.
which is the temperature distribution reached by the system without
heating effect. The system is then heated with a step input of the heat
generation to get an abrupt heating effect while the system is still
cooled by the coolant flow. After reaching a steady state. the heat

generation is suddenly cut off to get an abrupt cooling effect.

Figure 52 shows the temperature history of T; while the system
undergoes such a thermal process. and also indicates the times at which

the maximum heating rate and the maximum cooling rate are observed.

For each design change. three different values will be used in the
computer models. One is the original design. one is greater than the
original value numerically. and the remaining one is smaller than the
original value. Thus. with the corresponding computer-predicted
results. the effect of each design change on three design criteria (AT.
maximum heating rate and maximum cooling rate) then can be inspected.
From this information. a recommendation for an improved design of the

cryomicroscOpe conduction stage will be suggested.

104

onB<Hm<> mmemzdmdu mu 02H9<UHEMM>ZH mou Apv—B mm umDOHm

onedmmzmw B<mm who 9:0

mb¢m wzhgooo ZDZHx¢z

.mm.

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Sum. 93sz85
one¢mmzmu 9<am
as .3 .2: .8 .3 .3 4m \
d. ...N-
.
“ .37
7. .
. _
/ . ..2-
.
. . .
é. \ N, .3 a 333
(Ill-...“ \ /
me<m wzHadmm zzsz¢z
.9
.3
wine m> 1 mm:s(mmsxme

105
During the investigation of parameter variations. all the other
parameters are kept constant to assure the resulting effect is due to

the parameter of interest in the design change.

Tables 6 and 7 contain the numerical form of the results of all the
seven investigations. The graphs representing the trends of the effects
extracted from these Tables are shown in Figures 53-59 using a normal-
ized scale based on the original computer model. Specifically. the
normalized temperature difference. AT. will be defined as

AT = AT/AT (5.1)

original
where AToriginal is the steady-state temperature difference
between T, and T, of the original design. and
AT is the steady-state temperature difference resulting from

a smaller parameter value or a larger parameter value.

The reason for using the steady-state temperature difference is
because the steady-state temperature difference is the maximum tempera-
ture difference predicted by the model during the entire response for
all the parameter variation test conditions. However. this might not be
true for those test conditions where the system is initially at room
temperature. In these cases. there will be a time delay for the cooling
effect to reach the viewing hole region. and this time delay will cause
the maximum temperature difference to occur at some intermediate time
during the period of response instead of at the steady-state condiction.

Nevertheless. those test conditions for the parameter variations des-

Width of Viewing

Hole:

Thickness of
Glass Heater:

Magnitude of
Heat Generation:

‘1.3728 Watt

Width of Heat
Generation region:

0.15"
0.3"

0.6"

0.011"
0.022"

0.044"

0.6864 Watt

2.7456 Watt

0.73"

1.0"

1.3"

6.25
28.125

100.

38.
28.125

19.

20.0
28.125

44.35

33.04
28.125

24.44

Steady-State

Tymperature
T,-T, (C)

T,=—29’c
T,-3a‘c

T1=159°c

T,=62°C

T,=38°C

106

TABLE 6

Differencg

T,=14.3‘c

T,--so°c
T,=38‘c

T,=211°c

T1=9o°c

I;=3s’c

O
T1=o c

Maximum

Heating Rgtg
(

C/MIN.)

179.4
240.6

278.4

304.2
240.6

163.2

101.4
240.6

495.0

334.2
240.6

186.0

RESULTS OF INVESTIGATION OF PARAMETER VARIATIONS

"’5‘

x
O

( C/MIN

-255.0
-255.0

-246.6

-303.6
-255.0

-20508

-127e8
-25500

~510.0

-349.8
-255.0

-193.8

imum
in

Rgte
.)

107

TABLE 6 (CONT'D)

RESULTS OF INVESTIGATION OF PARAMETER.VARIATIONS

 

Steady-State Maximum Maximum
T}mperatur;)bifference Heaté7gIfigtg_ CQ?IC7SI§?;e
Thickness of
Copper Plate:
0.1122" 27.43 233.4 -251.4
‘o.224s" 23.125 240.6 -255.0
0.449" 26.09 240.6 -255.0
Material of
Metal Plate:
Brass 25.4 241.2 -252.0
Aluminum 27.27 241.2 -254.4
.Copper 28.125 240.6 -2ss.o
..Position of
Coolant Channel:
0.275" 28.69 240.6 -255.0
0.875" 31.03 240.6 -256.2
‘1.375" 28.125 240.6 -2ss.o

.Indicates original configuration

..Distance from Center of Viewing Hole to Center of Coolant Channel

THERMAL PROPERTIES OF MATERIALS FOR METAL PLATE

Thermal Conductivity Thermal Diffusivity
(W/M-°K) (MI/Sec)
Brass 81.36 2.48-10"
Aluminum 235.42 1.04.10“

Copper 484.4 1.32-10“

109

cribed in Sec. 5.1 are close to the actual situation. That is. before
ultilizing the cryomicroscope conduction stage to study the heat
transfer characteristics of biological samples. the system is cooled by
the coolant flow and reaches the steady-state temperature distribution.
The heat generation is then applied to the stage to control the tempera-
tures around the viewing hole region where the biological samples is
placed. So. the steady-state temperature differences are eligible to be
the observant quantities in determining the effect of parameter varia-
tions for the test conditions described in Sec. 5.1 (or Figure 52).

The normalized maximum heating rate and maximum cooling rate are defined

Heating Rate
Normalized Heating Rate = (5.2)
(Heating Rate)

 

original

and

Cooling Rate
Normalized Cooling Rate = (5.3)
(Cooling Rate)

 

original

where the heating rate and cooling rate actually mean the maximum
heating rate and maximum cooling rate for a particular
temperature history (see Figure 52). and the subscript
”original” indicates the corresponding values of the

original result.

The parameter value of each design change is also expressed in nor—
malized scale based on original design. The graph representing the
effect of each design change is explained below and will be discussed in

next section according to the Table 8 which summarizes the effects of

110

parameter variations.

Figure 53 shows the effect of the variation of the width of the
viewing hole. The original width of the viewing hole was 0.3" (0.15" on
the computer model diagram. Figure 6. which represents only the left
half portion of the symmetric cross-sectional area). The two additional
values used to observe the effect of this design change are 0.15". and

0.6" respectively.

Figure 54 shows the effect of the variation of the thickness of the
heating glass. The original thickness of the glass heater is 0.022".
The two variations of the thckness of glass heater are 0.011" and 0.044"

respectively.

Figure 55 shows the effect of the variation of the magnitude of the
heat generation. The original power input is 1.3728 Watt. Two varia-

tions for this design change are 0.6848 Watt and 2.7456 Watt.

Figure 56 shows the effect of the variation of the width of the
heat generation region. In this investigation. the total power input is
kept constant. and the width of the heat generation area for this con-
stant power input is changed. The original width of heat generation
region is 1" (0.5" in the computer model). The two variatons are 0.73”

and 1.3" respectively.

Figure 57 shows the effect of the variation of the thickness of the

copper plate. The original thickness of the copper plate is 0.2245"

 

 

STEADY-STATE

3. TEMPERATURE DIFFERENCE:
MAXIE‘IUI-I HEATING RATE:
MAXIT‘IUE COOLING RATE

2.

1.

O, +

 

 

O 1. 2.
WIDTH 0F VIEWING HOLE (NORMALIZED)

FIGURE 53 EFFECT OF CHANGING WIDTH 0F VIEWING HOLE

112

STEADY-STATE
TEMPERATURE DIFFERENCE:

 

MAXIHUM HEATING RATE:

2.
MAXIHUN CCOL HG RATE:

1.

 

 

 

0.
1. 2.

THICKNESS OF GLASS HEATER (NORMALIZED)

FIGURE 54. EFFECT OF CHANGING THICKNESS OF HEATING GLASS

113

STEADY-STATE
TED-IF ERATURE DIFFERENCE:

 

run-2m: HEATING RATE:
2. f9 "" '—
/ ammo cocune RATE:

1.

 

 

0. fig
0. 1O 2.

MAGNITUDE OF HEAT GE'Tnv‘ERATION (NORMALIZED)

FIGURE 55. EFFECT OF CHANGING MAGNITUDE OF HEAT GENERATIO‘;

11L:

 

 

 

 

STEADY-STATE
TEI'EPERATURE DIFFERENCE:
MAXIMUM HEATING RATE:

2. — — _ —
MAXIMUM COOLING RATE:

1.

O. #:—

 

005 1. 1.5
WIDTH OF HEAT GEI‘IEPATIOI-I REC-10:: (near-31.12313)

FIGURE 56. EFFECT OF CHAI‘lGIIJG WIDTH OF HEAT GENERATION REGION

115

STEADY-STATE
TEMPERATURE DIFFEREI‘JCE:

 

I'LAPZII-IUI'I HEATIIEG RATE:

1. -----------a "' — "
IRXIIJUI-l COOLIITG RATE:

0.5

 

 

Os 1. 9
THICKNESS OF COPPER PLATE (KORE-1.12.1221)

FIGURE 57. EFFECT OF CHANGING THICKNESS OF COPPER PLATE

 

116

STEADY-STATE
/TEMPERATURE DIFFERENCE

 

 

 

 

1° 9/3 43
Co I :
O. 6.16): COQEJCJ 1.
BRASS Al COPPER
THERI-IAL CONDUCTIVITY (IIORI-LILIZED)
I
Ia '..5.}".I .UI'. CCOL II G RATE
1. F — — — ..— -- b fi<3
'rJ—JI uUI; I-EATII..: A RATE
0. fi—
0. (3.135. 0.75., I.
BRASS Al COPPER

THERMAL DIFFUSIVITY (IIORW.I.—IZED)

FIGURE 53 EFFECT OF CHANGING MATERIAL OF METAL PLATE

117

STEADY-STATE
TEMPERATURE DIFFERENCE:

lflAXIZfUI-‘E HEATING RATE:
1. ...—"*1“... — _ — _

IIAXII-IUE'. COOLING RATE:

 

 

 

O. A '—
0. CQC 30;; Io

DISTANCE FROM CENTER OF VIEWING HOLE TO CENTER OF COOLANT CHANNEL
(NORMALIZED)

FIGURE 59. EFFECT OF CHANGING POSITION OF COCLANT CHANNEL

118

while the two variations are 0.1122" and 0.449" respectively.

Figure 58 shows the effect on the three design criteria due to
changing the material of the metal plate. Aluminum and Brass are used
here. because they are cheap and readily available. Since the compari-
son of temperature difference is based on steady-state temperature
difference. the thermal conductivities of these metals will be consi-
dered as the design parameter. On the other hand. the comparison of
heating rate and cooling rate is based on the dynamic characteristics of
the designs. so the thermal diffusivities of these metals will be used

as the design parameter in this case.

Figure 59 shows the effect of the variation of the position of the
coolant channel. The distance from center of the viewing hole to the
center of the coolant channel will be used as the design parameter to

indicates the position of coolant channel.

5.2 DISCUSSION AND CONCLUSION

Table 8 summarize the influence of each design change on the design
criteria. The sense of the vector indicates whether the influence is
positive or negative while varying the particular parameter.
Specifically. the vector which points upward at the first place of the
first column means that the temperature difference. AT. will increase
while increasing the width of the viewing hole. The notation "--" in
Table 8 means that there is no significant influence on the design cri-

teria by varying the corresponding design parameter.

THE SUMMARY OF PARAMETER VARIATIONS

Steady-State

Tgmpgrgture Qifference

Width of Viewing
Hole

Thickness of
Glass Heater

Magnitude of
Heat Generation

Width of Heat
Generation region

Thickness of
Copper Plate

Material of
Metal Plate

Positon of
Coolant Channel

?

Maximum
Heating Rgte

f

Maximum
Cooling Rgte

120

From row one in Table 8. it may be concluded that AT and the
maximum heating rate will increase. and the maximum cooling rate will
not be affected while changing the width of the viewing hole. To see
this. refer to the computer model diagram (see Figure 6). If the width
of the viewing hole increases. then the distance between T3 and Ta will
increase. because T1 is selected as the temperature of the center of the
viewing hole, and Ta is selected as the temperature on the edge of the
viewing hole. This increment in distance results in a larger tempera-
ture difference between T1 and Ta while keeping all the other parameters

constant.

As mentioned in Section 4.2.3. the contact resistance between the
glass heater and the Kapton film is smaller than the thermal resistance
of the air. and the thermal energy generated inside the viewing hole
area is difficult to diffuse. So. a large viewing hole also means the
thermal energy generated tends to accumulate more. and this also results
in a larger temperature difference between T1, the temperature of the
center of viewing hole. and T3. the edge of the viewing hole. By the
same reason. while thermal energy is generated the temperature inside
the viewing hole will increase rapidly because only a small amount of

thermal energy will be removed per unit time.

From another point of view, the variation of the width of the view-
ing hole changes the local geometry around the location of interest. TH.
Consequently. this variation has a significant influence on the local
heat transfer characteristics. Since the maximum cooling rate occurs

after the heat generation has been cut off. the thermal energy rapidly

121
diffuses from the viewing hole to the coolant channel. The possible
effect on this cooling process by changing the width of the viewing hole
is that it determines the area which is still in contact with the Kapton

film so that can be used for passage of this energy diffusion process.

It is possible to imagine that if the variation of the width of the
viewing hole is not significant compared with the total area which con-
tacts with the Kapton film. then the influence of this change on the

maximum cooling rate will be trivial.

Furthermore. Figure 53 reveals that the effect of this design
change on AT is significantly greater than on the maximum heating rate
and the maximum cooling rate. It means that by deceasing the width of
the viewing hole. the horizontal uniformity of the temperature distribu-
tion within the viewing hole region can be significantly improved
without losing too much ability in conducting the experiment which
requires rapid responding system. So. this will be a useful information

for optimal design consideration.

From the second row of Table 8. it is observed that AT and the max-
imum heating rate and the maximum cooling rate all will decrease while
increasing the thickness of the heating glass. Generally, the thermal
energy generated will diffuse to both sides of the heat generation plane
(see Figure 6). The amount of thermal energy diffused to either side is
determined by the relative thermal resistance of the materials in either
side. Figure 60 shows a descriptive diagram of this thermal resistance

for a simple one-dimensional case. If [left is smaller than Kright'

122

Qleft + Qright = HEAT GENERATION

Qleft Qright

‘—

"‘"""\/A\v/\vr :' “v/A\v/A\V/——_—*
.l. .1.
Kleft Kright

 

HEAT GENERATION

 

Kleft Kright

 

 

 

 

FIGURE 60. THERMAL RESISTANCE DIAGRAM FOR A COMPOSITE
MATERIAL WITH HEAT GENERATION

123

then most of the thermal energy will diffuse to the left hand side.

Applying this model to the case of the cryomicroscope conduction
stage. an increase in the thickness of the glass heater results in a
wider passage for the diffusion of thermal energy to the coolant chan-
nel. This. in turn. means the thermal resistance in the lower side of
the heat generation plane is decreased (see Figure 6). Thus. the ther-
mal energy diffused to the upper side of the heat generation plane is
decreased. and this tends to flatten the temperature distribution on the
top of the cover glass which is on the upper side of the heat generation
plane. This also results in a decrease of the temperature difference

between T1 and T,.

Since the heating rate and cooling rate are dynamic characteristics
of the system. it will be better to explain the behavior of the heating

rate and the cooling rate in terms of thermal diffusivity.

Increasing the thickness of the glass heater always increases the
thermal capacitance of the glass heater. With a larger thermal capaci-
tance. it tends to slow the heat transfer response to a thermal
stimulation. It also means the heating rate and cooling rate will slow

down.

From the third row of Table 8. it is apparent that AT and the maxi-
mum heating rate and the maximum cooling rate will increase while
increasing the power input. It is easy to see in this case that AT and

the maximum heating rate will increase for a larger power input as com-

121+

pared with a smaller power input. As to the cooling rate. since the
system has a larger power input. it also means the steady-state tempera-
ture will be greater than that of a smaller power input. Thus. there
will be a larger temperature gradient after cutting off the heat genera-
tion. and this causes a faster cooling rate as compared with a smalle

power input.

From the fourth row in Table 8. it can be seen that AT. the maximum
heating rate. and maximum cooling rate. will decrease while decreasing
the width of the heat generation region and keeping the power input con—
stant. Since power input is kept constant. to widen the width of the
heat generation region is the same as decreasing the heat generation per

unit area.

From the local point of view. this is equivalent to decreasing the
power input. Consequently. AT. the maximum heating rate. and the maxi-
mmm cooling rate will decrease. This reveals that the heat transfer
characteristics of T; is affected by the local heat generation rather

than the heat generation area.

The fifth. sixth, and seventh rows of Table 8 indicate that AT. and
the maximum heating rate. and the maximum cooling rate are not signifi-
cantly affected by the geometry. material of the metal plate. and the
coolant position on the metal plate. This can be account for by the
high thermal conductivities of metals (Copper. Aluminum. and Brass).
The metal plates (Aluminum. Copper. and Brass) all function with approx-

imately the same efficiency in transporting the accumulated thermal

125

energy within the viewing hole to the coolant channel.

Figures 61 and 62 show some comparisons of the results of the
present investigation of the parameter variations with the results of
Khompis's work (Tables 7.1 and 7.2 of reference [3]). The normalization
of the coordinate axes on these figures are based on Ihompis's defini-
tion of normalization. It appears that there are significant deviations
between the results of the present research and Khompis's results. This
is because the cryomicroscope conduction stages used and the computer
models for these two independent research efforts are basically dif-

ferent.

Furthermore. the curve of the present research in Figure 61 has a
slightly positive s10pe. This means AT. the horizontal temperature
difference between T1 and T,, will increase as decreasing the width of
the heat generation region and this is contradictory to the other curves
in Figures 61 and 63. This is because of the improper method of the
normalization. However. it is still possible to get a meaningful com-
parison of these two results. if an appropriate method of normalization
is employed. Figures 63 and 64 show a comparison of the results of the
investigation of parameter variations with Khompis's results with a nor-
malization which is based on the definition described in Sec. 5.1.
Figures 63 and 64 reveal that the results for the present case and the
previous case are quite consistent for the variation of the width of the
heat generation region and the variation of the thickness of the glass

heater.

126

RESULT OF PRESENT

 

 

 

 

RESEARCH:
20.
KHOHPIS'S RESULT:
15.
2.? f
0
g 10. I
m 1
1 A\
l ‘1
.~‘~
50 \A‘
“A
0' :
O. I. 2.
WIDTE OF HEAT GEIERATION REGION (NORMALIZED)
T - T
men (as) =—‘———i
T1 " TC

FIGUREEH. COMPARISON WITH KHOKPIS'S RESULT OF CHANGING THE
WIDTH OF HEAT GENERATION REGION (KHOMPIS'S DEFINITION
FOR NORMALIZATION)

W

1277

RESULT OF PRESENT

2. RESEARCH:

 

KHOMPIS'S RESULT:

)

ERROR (?

 

 

1Q.
\
‘RTA-.
50 \1\\
‘15
O. -*
O. 1. 2.

TRICKNEss 0E GLASS HEATER (NORMALIZED)

FIGURE: 62. COMPARISON WITH KHOHPIS'S RESULT OF CHANGING THE
THICKNESS OF GLASS HEATER (KHOMPIS'S DEFINITION
FOR NORMALIZATION)

STEADY-STATE TEMPERATURE DIFFERENCE (NORMALIZED)

2.

 

128

RESULT OF PRESENT
RESEARCH:

 

KHOHPIS'S RESULT:

1. 1.5
WIDTH OF HEAT GENERATION REGION (NORMALIZED)

FIGURE 63 COMPARISON WITH KHOMPIS'S RESULT OF CHANGING THE
WIDTH OF HEAT GENERATION REGION

TEMPERATURE DIFFERENCE (NORMALIZED)

STEADY-STATE

129

RESULT OF PRESENT
RESEARCH:

 

KHOHPIS'S RESULT:

2.

1.

 

 

O. 1. a.
THICKNESS OF ;LLJS HEATER (NORMALIZED)

FIGURE 6#. COMPARISON WITH KHOMPIS'S RESULT OF CHANGING THE
THICKNESS OF GLASS HEATER

130

As mentioned in Chapter 1. an ideal cryomicroscope conduction stage
should possess an uniform temperature distribution within the viewing
hole and be capable of a fast thermal response. This indicates the need
for three design criteria. In particular. it is desireable to make
design changes which minimize AT while maximizing the maximum heating
and cooling rates. A close observation of Table 8 shows that all of the
variations of parameters lead to contradictory influences on the three
design criteria. That is, minimizing the temperature difference will

also bring about a minimal and not maximal thermal response.

A compromise has to be made in order to decide which criteria has
priority. In this case. the uniformity of the temperature distribution
within the viewing hole area is considered to be most important. Based
on this criteria. a recommendation for a better design of the cryomicro-
scOpe conduction stage can be made according to Table 8. The
recommended improvements include:

. Small Viewing Bole
. Thick Glass Heater
. Small Heat Generation

. Wide Beat Generation Region

Reviewing Figures 53 to 59. each design change has similar quanti-
tative (normalized) effect on AT. maximum heating rate. and maximum
cooling rate except the case of changing the width of the viewing hole.
The increase of the width of the viewing hole leads to an increase by
square of AT while the maximum heating rate and the maximum cooling rate

approximately keep constant. This leads to an important conclusion that

131
reducing of the width of the viewing hole is the major consideration for
an improved design for the cryomicroscope conduction stage. The materi-
al and the thickness of the metal plate are trivial effects. and the
coolant channel can be located at any positon on the metal plate.

Figure 65 shows a descriptive scheme of this recommandation.

5.3 SUGGESTIONS FOR FUTURE.STUDY

Several suggestions for future study are stated below.
A. Investigate the temperature distribution in the Y‘direction (see
Figure 3). which includes the important information about the error
between the temperature of temperature sensor. Ti, and temperature of
biOIOSiCSI Samples. 1;. This can also be used to find out the lag times
between the heat generation and the effect at the temperature sensor
which is an important consideration for control stability.
B. Investigate the effect of changing the thickness of the cover
glasses on the temperature difference between T1 and T;. C. Conduct
experiments to varify the results of parameter variations made here.
D. Investigate the effect of changing the magnitude of the heat
transfer coefficient of the coolant flow and the effect of the contact
resistance.
E. Simplify the computer model and examine the possibility of eliminat-
ing some unnecessary considerations. such as free convection on the
vertical boundary because of the small dimension in the vertical direc-
tion. or consider a lumped model for the metal plate because of the
essentially uniform temperature distribution within it.

F. More work on pre-processing and post-processing software develop-

132

ZOHB<szzzoomm MES mo mzmmom m>HBmHmommn omm mmeHh

MACE GZHBHH> Addzm

 

_ me¢qm A<Bmz sze

 

 

 

 

 

 

_ qmzz¢mu _
_ ezeqooo _

.r////// ona¢mmzmu a<mm

AA<Em 38H; onwmm
ZOHB¢mmzmo B<mm mQHB

mme¢mm mw¢dw MUHmB

 

 

 

 

 

 

 

ment.
G. More work on the Textronic 4114 type of display and/or color coding

of temperatures.

APPENDIX A:

CALCULATION OF PARAMETER VALUES

This section.will describe in detail how the parameter values used
in the computer model were evaluated (see Tables 3 and 4). Some of
these parameters are already documented. some have to be calculated .
and the rest of the parameters are estimated values assigned from rea-

sonable estimates.

Calculation of the parameters of case 4 will be used to illustrate
the process of finding these parameter values. The process of finding
the parameter values for other cases will be similar to the process des-

cribed below.

A. no and ac: The thermal conductivity and diffusivity of cOpper.

They are documented and can be found in reference [5].

B. Is and as: The thermal conductivity and diffusivity of glass.
The material of the rectangular glass heater of the system is Borosili-
cate glass. and it is a product of Corning Glass Works Company . and
officially named 'Electrically Conductive Coated Glass'. The thermal
properties of Borosilicate glass can be found in reference [5]. To sim-

plify the computer model. the same property values will be used for the

cover glasses used in this system.

134

135
C. Kf and of: The thermal conductivity and diffusivity of Kapton
Polyimide film. which is a product of Du Pont Company. The thermal pro-
perty information offered by the manufacturer are available only at
certain temperatures. Thus. extrapolation of these data has been done
by the investigator to get the thermal property values for the tempera-
ture range of the experiment. Figure 66 shows the extrapolation for

thermal conductivity of Kapton film.

D. K. and “e: The thermal conductivity and diffusivity of glue
epoxy. The glue epoxy is used to bond the glass heater and cover
glasses together. Since the thermal property information is not avail-
able. the property values of Teflon (Polytetraflouroethylene. which can
also be found in reference [5]) are used as estimated values for input
to the computer. By introducing a correctional factor. the proper
values which make the numerical results match experimental results then

can be found by trial and error on interactive computer terminal.

E. Xi and “i: The thermal properties for the imaginary material
which models the contact resistance between the glass heater and the
Kapton film (see Figure 3). Since this layer of material does not actu-
ally exist. the property values can be chosen freely as long as they can
make numerical results match experimental results. The proper values

are found on the interactive computer terminal by trial and error.

F. G: The normalized heat generation per unit area. In case 4

(see Table 1). the power input is

p = v’ / R = (7.66 volt)’ / 54 Ohm = 1.0866 Watt (A.1)

THERMAL
CONDUCT
(W/Ho'K)

0.19

0.18

0.17

0.16

0.15

0.13

 

136

TEMPERATURE(°K)

 

T

100 200 300 400 500 600

FIGURE 66. EXTRAPOLATION FOR THERMAL CONDUCTIVITY
OF KAPTON FILM

157
where P is the power input. and V is the voltage across the glass
heater. R is the electrical resistance across the glass heater. and G

is defined as

P L
G = _ o—i (A02)
A KET,

In this equation. L: is the basis for the normalization of the X coordi-
nate. and it is selected as the total width of the computer model in the
present study (1.5". see Figure 6). Is is the thermal conductivity of
the glass heater. and T. is the basis for normalization of the tempera-
ture scale. Since the operating temperature range of case 4 is from
-l96°C to 26°C, which 1. approximately 250°C. 250°C is used for the
basis of normalization for the temperature scale. The group KgTOILx
then serves as the basis of the normalization of the heat generation per
unit area. In the above definition. A is the effective area of heat
generation which excludes the region covered by silver epoxy. besause
the silver epoxy is such a good conductive material that the electric
current will pass through the silver epoxy instead of the metallic oxide
coating. Then 6 becomes
1.0866W 1.5"

. g 0.27
(1")’ 250°K-0.95W/M-°K

 

 

6- Ht: This is the heat transfer coefficient for free convection
on the upper surface of the cover glass (see Table 1). According to
Sec. 7.3 and 7.4 of reference [6]. the empirical relation for the Nus-

selt number is

138

an = o.5s(crf-Prf)°-‘ (A.3)

where Grf is the Grashof number and Prf is the Prandlt number.
The product Grf-Prf is the Rayleigh number. Raf. and is evaluated at the
'film temperature'. T}. which is the arithmetic mean of the wall tem-

perature . T

V. and the ambient temperature. 1;.

Tf =(T¢+T;)/2 (A.4)
In case 4. the film temperature is
Tf = (296°x+155°x)/2 = 225°:
where the average temperature of the upper surface of the cover glass is

chosen as 155°R. According to reference [6]. the Rayleigh number of air

at 225°: is

R.f (Tg-T;)-x;-Prf-(g-p.p’/p’)

(141’s)-(o.875")'-o.734-(ss.1-10' °r“l”)

6.6.10‘

where Xe is the characteristic length and is chosen as the width of the
cover glass. 0.875". in this case. The values for the dimensionless
group g-fl-p’lp' and Prandlt number. Prf. can be found in Table Arl of

reference [7]. Substituting the calculated values into equation (A.3).

139

the average Nusselt number is

N3. = (o.ss)-(660110>°°’ = 9.203

From the definition of the Nusselt number. the heat transfer coeffi-

cient. Ht. can be expressed as
-—- a 0
Et = Nuf-K/Xc = 7.68 W/M - K (A.5)
This is the calculated value for Ht' and it can be modified by
introducing a correctional factor so that the numerical results can
match the experimental results. In case 4. the correctional factor
found by computer matching for at is 0.7.

H. Eb: This is the heat transfer coefficient for free convection
on the lower surface of the cOpper plate. An empirical relation for
this geometry also can be found in Sec. 7.3 and 7.4 of reference [6].
It is

Nuf = 0.15(Grf-Prf)°°” (A.6)
The film temperature is

Tf = (173°r+102°K)/2 = 133’s

The Rayleigh number. Raf. at 138°K is

1A0

(70°x)~(2.s")’-0.769-(s.7-10' °x“u")

Raf

12.10'

Substituting this value into equation (A.6). it becomes

Nuf = 0.15-(12-10’)'-" = 74

. and the heat transfer coefficient is

ab = 14.75 W/M’o'K

1. 3,: This is the heat transfer coefficient for free convection

on the vertical boundary of the copper plate. The film temperature is

Tf = (173°x+102°x)/2 = l38°K

. and the Rayleigh number is

Raf = (70‘x)-(0.0245")’-0.769-(s.7-10' °x"u”)

= 112.9

According to Figure 7-7 of reference [6]. the corresponding Nusselt

number is

141

. and the heat transfer coefficient is
H, = 54.3 wvu’-°:

1. HI: This is the heat transfer coefficient for free convection
on the vertical boundary of the cover glasses and epoxy. The film tem-
perature is

Tf = (269°:+143°:)/2 . 220°:

. and the Rayleigh number is

(153°K)'(0.014")'-O.735-(63.66-10' °:"u’°)

Raf

= 3.22

According to Figure 7-7 of reference [6]. the corresponding Nusselt

number is

. and the heat transfer coefficient is

H1 = 98.35 Wlu’-°K

K. Be: This is the heat transfer coefficient for forced convec-

tion due to coolant flow. Figure 67 shows the sketch of the path of the

142
coolant flow. Tb calculate Bc' the coolant flowrate has to be calculat-
ed first and substituted into the equation of Nusselt number to
calculate the required heat transfer coefficient. Be. There are several
assumptions which should be mentioned before starting the calculation.
1. The bending effect of the flow path is neglected.
2. Bernoulli equation is assumed to be applicable from planes (J)
(I) to on the sketch.
3. The effect on effective orifice area of fluffy ice gathering
around the orifice of the plenum due to low coolant

temperature inside the plenum is neglected.

4. The flow velocity inside the plenum is neglected. that is.

5. The coolant is assumed to be incompressible only between
planes (1) and (R). that is. pj 2 pk.

6. Assume no leakage occurs along the flow process. that is.
all the mass flowrates are equal. M1 = Mj = Mk.

From Figure 67. the Bernoulli equation for flow between planes (1) and

(K) is
p v? P V’
k k
91 2 pk 2

Rewrite the equation. it becomes

 

3
Vt = Pj-Pk g Pj-P‘tm a AP
2 "1: P]: P:

Vk = ZAP/pk (A.8)

1A3

301E 529000 ho FEE .50 mmeHh

5.5.6535 5.5.3.;

x mzfim

w mzmmm

 

mUHmHmo _

 

zmwomBHz 20mm

 

23.5an

 

_
1L
.
_
_ mz<e
_
_
_
_

 

 

 

 

H mz<qm 11;..1 1: x..: I. I..J.n%.x

rv rill/Izmmzzemu

8733000

ace mmpmmmmm
$£>9fi£

 

 

 

 

 

 

 

Egan ZOHBUDQZOO
mmoomomonowmo

w.

The mass flowrate is

Mk = 0.61.51kak (A.9)

where 0.61 is the flow coefficient for square edged orifice (reference

[8]).

From the experimental data. the difference in water column heights

of the U-Thbe is 1.5". that is.

AP = 0.0158 psi.

From Table A.5.2 of reference [9]. the density of Nitrogen at 143°K
is 0.15 lbm/ft'. Substituting these values into equation (A.9). the

mass flowrate is

“h = 6.15 lbm/hr. = M1

An average velocity. V. is defined as

1'

Mi 3 PiAivi (A.10)

where Ai is the cross-sectional area of the coolant channel. and pi is
the density of two-phase flow inside the coolant channel. At this
point. a problem is encountered. That is. the 'quality' of this

two-phase flow. Xq. is not available. The direct way to eliminate this

145

problem is to quess a reasonable value for it.

Since the quality. Xq. is apparently greater than 0.5 and also is
certainly smaller than 1 by definition. the quality is selected as mid-

way of these two values. that is.

X = 0.75

From Table A.5 of reference [9]. the density of Nitrogen at 14.7 psi.
with quality 0.75 is 0.3826 lbm/ft'. then the average velocity. 71, will

be

‘71 = (0.3826 lbm/ft')°(0.0205 1n°)/6.15 lbm/hr. = 31.44 ft/sec.

Raving calculated the required average velocity. the next step is to
calculate the heat transfer coefficient for forced convection of
two-phase coolant flow at -196°K. To simplify the problem. the Nitrogen
will be treated as at a saturated gaseous state in the following pro-

cedure. The Reynold number. Re. of this flow is

Re = 4-Xc-V-p/u = 31820.6 > 2300 (A.12)

where Re is the characteristic length of the channel. It appears that

the coolant flow is turbulent flow.

Since the Reynold number is greater than 2300 and the Prandlt

number is less than 1. the empirical equation (13.9) of reference [7]

146

can be used. It is

Nu = 0.022(P:)°-‘(Ro)°-' (A.13)

= 0.022(0.787)°-‘(31820.6)°-' = 78.07

. then the heat transfer coefficient is

B0 = 182.2 WIu°-°x

This value is taken as a reference value so that it can be input to
the computer for finding the proper one. The resulting value for RC is

583. W/M’o'K. and the correctional factor is 3.2.

L. Ebo: This is the heat transfer coefficient for the lower sur-
face of the heating glass which is exposed to the air through the
viewing hole. The reason for distinguishing this parameter with that
for the lower surface of the copper plate. Eb. is that they have dif-
ferent geometries. Furthermore. the heat transfer coefficient below the
viewing hole. Hbo’ is expected to have a significant influence on the
temperature distribution around the viewing hole which is of concern in
the present study. Basically. Rb will be used as a reference value for
abo' and the more prOper value will be found on the interactive terminal
by trial and error. In cases 3.4 and 5. no correctional factor is need-

ed for Hbo‘ but in cases 1 and 2 the correctional factor is 0.5.

11:7
M. Hag‘ This is the heat transfer coefficient used to take care
of the effect due to the contact of the metal conductor with the silver
epoxy (see Figure 2). Thus. the value can be chosen freely so that the
numerical results can match the experimental results. at is used as the

reference value and modified with correctional factor to get a proper

value. In case 4. the correctional factor is 4.

 

APPENDIX B:

DERIVATION 0F FORWARD-DIFFERENCE EQUATION FOR INTERIOR NODE

An interior node. which is located neither on the boundary nor at
the corner of a forward-difference griding diagram. is shown in Figure
68. According to chapter 8 of reference [10]. the energy balance for

this rectangular region centered by the node (i.j) is

on + an - on - Q“ + 0y1 - 0y2 + g .. H pcaT/atdv (13.1)

where Q's are heat fluxes. and g is the total heat generation inside the
nodal system. The C inside the integral is the specific heat of the
material. and dV here is the differential volume. so the integral stands
for the rate of change of the internal energy of the nodal system.

These heat fluxes can be approximated by

 

 

 

 

 

 

AT T“ . - T‘L .

Q11 = "Kg—f' ' 1“, i 1,J (3.2)
2 A11
AY Tn — Tn_ .

sz = -K1 a. 1.1 i 1,] (8.3)
2 AX;

9!. = -K3_-: 0 + ,j i,j (3.4)
2 8:,
AI '1“ - Tn .

ox. = -:1-—:- “1'1 i'J (B.5)
2 AX,
AX: + Ax: Tgpj — 1?.j-1

0.]. = -:1 - (8.6)

m9

 

'-—-dX1 i v—-dX2

 

)

(I-I,J+1) (I,J+1) (I+1,J+1)

. dYZ {'7/' I Q K2,C2,P2
QX] 1 x3
‘ ' \q
(I-1,J) N (I+1,J)
.__*J\— |-—-—-1—-
QXZ L— — Qxh

\‘y
A

 

dY1 K],C1,P1
Qy]
I O
(I-1,J-1) (I,J-1) (I+1,J-1)

FIGURE 68. INTERIOR NODAL SYSTEM WITH ENERGY TERMS INDICATED

150

n
AX: + Ax: T1.3+1 ’ TI.j

 

03,, = -x, (B.7)
2 AY

Assuming that the materials are isotropic so that the specific heat
C and the density are constants in all directions. and assuming that the
rectangular region is small enough so that it can be treated as a lumped
system. the rate of change of internal energy is

AY1 Al: + AX. TTT} - 1?,3
ijca'r/atdv = p1C1 . .

 

 

2 2 At
AY AX + AX Tn+1 - I“
a a 2 1,1 1,1
+ 03C, °:~ (3.8)
2 At
By defining

RX = szle1 (B.9)
RY = AYzlAY1 (3.10)
R: = Ira/x1 (3.11)
R“ = azlal (3.12)

. and substituting all of these_equations into equation (8.1). the equa-

tion becomes

. n+1 -
(1+2: RK/Ro.) TM. T11 23

20, At AX,AY3K,(1+RX)

 

 

 

(1+RK-RY) Til-1,3 - (1+1/errj1‘,j + (Illuirr'i’fl,j

 

 

 

1+Rx (AX1)3

151

In - (1+RK/RY)T“ . + (ax/mm“
is -1 a a
+ 3 i J i 5+1 (8.13)
(AT1)’

 

To normalize the equation. the normalization of several variables are

defined as
Ai’= AX/Lx (B.14)
A? = AY/Ly (B.15)
EL = Lx/Ly (8.16)
E = s/ToK1 (B.l7)
A? = At.61/(Lx)° (8.18)

where (Lx)a/u1 is the time constant of material 1.

In the definition of AT. the time constant of material 1 is used as
the basis of normalization. When there are several materials composited
together. this definition will raise some problems. because the bases of
normalizations will be different for different layers of composite

materials. Thus. Rt is introduced and is defined as

Rt 3 teui/(Lx)z (3.19)

where to is a selected global base for normalization of AI.

152

Since the operating temperature range during the experiment is from
~196°C (coolant temperature) to 26°C (ambient temperature). approximate-

ly 250°K of span. the normalized temperature can be defined as
0 = (T+200°:)/'r. = (T+200°:)/250°: (8.20)

where T is in °C.
Applying these definitions. the forward-difference equation for an

interior node becomes

(1+RK-RY/Ra) eff} - 9111.3 2281.

 

 

2:: A? (1+RX)A-X'1AT1

+ (1+RK-RY) 03.1.1 - n+1/Email. + (unwind

1+RX (AX-1) ‘

 

 

cg, j-1 - (1+RX/RY)0‘£'j + (RR/RYWf’ 3.1

+ (RL)° (8.21)

 

(AYQ’
After slight modification. the equation then becomes the form appearing
in computer subroutine PLT.FDB. which contains all the thirteen nodal

equations.

APPENDIX C:

FLOWCHART‘S AND LISTING OF COMPUTER PRmRAM

To aid the explanation of the computer program. the flowcharts for
most of the computer program are shown in Figure 70-76. Immediately
following these flowcharts are the program listing. Brief comments are

included in the listing to aid the reading of the program.

The name of the main program is MBPLT' which selects the output
forms and assembles forward-difference equations for the system by cal-
ling individual forward-difference equations from PLT.FDE. which
contains thirteen subroutines corresponding to the thirteen nodal equa-
tions. Because of the simplicity of these subroutines. only the

subroutine for the interior node. INNER. will be presented.

Subroutine INOUT takes care of the input and output actions.
Subroutine GRAPBZD draws several charts which contain temperature histo-
ry information of T,, T,. T,. T}. and T,. and also displays the
temperature difference between T, and T, on a separate chart.
Subroutine GRAPBBD draws normalized three-dimensional temperature dis-

tribution profiles.

The logic of the group of subroutines ISOTHERM. SIDET. SIDMB,
SIDEL. SIDER is somewhat complicated. so it will be more clear to use
some verbal descriptions instead of presenting the flowchart only. The
purpose of these subroutines is to draw isotherms on the cross-sectional

configuration of the cryomicroscope conduction stage to indicate the

153

151+

temperature gradients within this system.

First of all. the irregular shape of the configuration is divided
into several rectangular regions. Then. for each rectangular region. it

is further divided into small grids such as shown in Figure 69.

The pragram will go around the boundary of each region to find out
whether the specific temperature of a isotherm is on it or not. If it
did find a point on the boundary which is at the specific temperature of
this isotherm. the cursor will move to this point. Starting from the
small grid containing this point on the boundary. the program draws the
isotherm by subsequently passing the cursor through small grids.
Usually. when an isotherm passes through a small grid. there should be
two intersections with the boundary of this small grid. So. the program
continues to find the other point on the boundary of this small grid.
After the cursor connecting these two points. the program proceeds to do
the same thing for its neighboring grid. Thus. the program connects the

isotherm with piecewise linear segments which pass through small grids.

Until the isotherm terminates at the boundary of a region. this
procedure continues for a different isotherm or a different region.
With the aid of four subroutines SIDET. SIDEB. SIDEL. and SIDER.
subroutine ISOTRERM can accomplish this task. The function of each

subroutine is stated below.

Subroutine SIDET is used to find the 'exit' point of the isotherm

segment which passes through a small grid and is entered from the upper

GZHBEOAQ :mmmeomH mom monmH> Q .00 mmmeh

 

 

 

 

 

L
w
W

155
‘LL
I
ll
I'll
III _
III
III

 

 

m
H /

M 1 /

onwmm

 

 

 

II

 

 

 

mzmmmBOmH

156
boundary of this grid. Subroutine SIDBB is used to find the 'exit'
point of the isotherm segment which passes through a small grid and is
entered from the lower boundary of this grid. Subroutine SIDEL is used
to find the 'exit' point of the isotherm segment which passes through a
small grid and is entered from the left boundary of this grid.
Subroutine SIDER is used to find the 'exit' point of the isotherm seg-
ment which passes through a small grid and is entered from the right
boundary of this grid. Because of the analogy in logic of these four

subroutine. only the flowchart of SIDRB will be presented.

The referencing relation of these program units is expressed with a

'tree diagram' as shown in Figure 77.

The subroutines of the standard FORTRAN library and the package
PLOTlO developed by Tektronix Company referenced by the present program

are omitted in the flowcharts.

157

 

 

 

 

 

 

".' V. r‘ 1*
10L d-U
‘."‘:“1.a
AiKJ‘AJiJ

CI'IVFF
'Iaf-‘e -J.(_ ./

INITIALIZE E-D
GRAPEIC ROME

SUI- . GPaAADZIED)

    
   

n ' BY:
-.J Wu).

 

 

 

 

TPUT?

    

 

f r1fi‘f“

. INI
13);. DIST
IIOQ, A133
UT ILIIIAL
P. DIST.

(EFIHRI' C17?

.24! ‘VOA

 

 

 

TIAL
REFU-

\Jer

T:II‘::
TIT? As] .)

 

 

 

 

 

IRITIALIZF
ITER
ABLES AND

ULITE EIOT

IETERS

 

COU-

,EIAC,VARI-

CALC-

1n

AI ’.

bIED OthR FARR-

 

 

 

 

 

 

SIART TIME
FfBHEHT

UZ‘ 1'

INC-

- 4.11

 

+ 3T

ICCUNTzICOUNT+1

 

 

<. ’V’.’
*1. rs

- *UL/“k—J

. C’uvfivu:aR'—'—'
. -JC’ elVa.&A

 

OF HAIR PROGRAM ‘MSPL T'

158

CALCULATE TEMP.
DIST. BY FORWA-
RD—DIFFERENCE ‘
METHOD (SUB. I-
NNER,SUP,SBTM,
SLEFT,...)

AND STORE IN
'UO(ILJ)W

 

 

 

 

  
 
 

CUT OFF HEAT

NO YES GENERATION

 

 

 
 
  

 
 
 

IT TIME TO OUT-
‘UT DATA?

   
 
 

T’UT DATA?

 

YES

   
   

OUTPUT NUMERIC-
0 III: TEI’IP. DIST.
AND RATE OF TE-
MP. CHANGE
(ENTRY OUTPTEM)

   
  
 
 

 

 

1

NO WANT YES

  
 

 

TPUT?

 
 
  

    

IS
IT TIME TO PLO
-D GRAPH?

   
 

IS
IT TIHE TO PLO
ISOTHERHS?

NO NO

 

           
   

 

 

 

 

 

 

 

 

 

 

 

 

 

‘ YES
max: ISTHERMS DRAW 3-D TEHP.
on 2-9 CRARH DIST. PROFILE
gggk. ISOTEE- (ENTRY GBDENTI)
.. ‘R -
EU; I J

FIGURE 70. (CONT'D.)

159

0

593 GRAPHIC OU- -
TPUT? 1

YES

 

 

ILITlALIZE 2-D
NO GRAPHIC MODE
(SUB. GRAPHBD)

’ I

 

\J‘J

WANT YES

HAFTS OUTPUT
?

NO

 

 

 

PLOT CRAFTS CF '
TEMP. HISTORIES

FERENCE

(EKTRY GZDENTZ)
STOP \—i l

 

 

 

 

FIGURE 70. (CONT'D.)

 

OPEI FILE UNITS
FOR I/O

 

 

 

READ IHPUT DATA
AHD WRITE THEM
ON OUTPUT FILE
FOR KEEPIKG RE-
CORDER

r l

 

:13:
6»;
.,<

 

 

RETURN

FIGURE '71.FLOWCHART OF SUBROUTIKE 'INOUT'

(71"... ".!f‘_"?“‘1
I51 1 VI. 1

I -.

-“"\"\ 1‘ ‘1
VI"’[’A\ F

L'-.vvi -. 1.1 .JH

.F’.T]FT1Y‘)F7‘ "I?"
\zL/ J .> J L;-‘ -‘-———-----

 

 

CEAEGE NORMALI—
ZED TEMP. AND
TIME IKTC REAL
SCALE

STORE REQUIRED
IhFORMATIOH FOR

TEMP. HISTORIES
DRAWING

OUTPUT NUMERIC-

 

 

 

 

 

 

 

 

 

CALCULATE RATE
OF TEMP. CHANGE
END OUTPUT

 

 

RETURN >

 

[CLOSE I/O FILES]

RETURN >

 

 

 

FIGURE 71.(CONT'D.)

 

* f': i"'-- T.“='*'N:‘ar *‘3’: T ”7}

V. -'_ ."md J ..uv 1.211 :L, J".‘------
*1" 1”“
H . .

I

INITIALIZE
GRAPHIC HOD

(:vRETURN:>

SECOID ENTRY PCIHT‘

 

—D

tr} \JJ

 

 

 

 

 

 

DRAW 5
F

Pflu CD
“A. a

b-fi
u-I
'31 .

:r» U
*u
m
*2:
m
h 1
I

 

 

 

DRAW 5-D TEMP.
DISTRIBUTION
PROFILE

 

 

 

{RENE-raj

 

 

CLOSE OUTPUT
CRAPEIC FILE

 

 

9:4
*—
C
(.
111
L4

\1
l\)

o
|>11

VEOWCHART OF SUBROUTINE 'GRAPEED'

TI'DNIT‘. Y—“ ff‘r‘Y"., “ "a :‘T’V‘.
' . ‘ v. , I I
‘ JAM}. L-"L‘\J _‘ \J-LA‘JI-------

 

 

 

 

IRITIALIZE 2-D
GRAPHIC MODE

 

 

 

 

(:RETURH >
EHTRY POINT

\.;2IJ ‘7

 

 

 

 

   

DJAW BASIC COH-
FIGURATIOX OF

.-7“ évcwr”
TEAL U IIUl ul'.

 

 

 

 

THIRD ENTRY POINT

 

 

SET UP THE WIN-
DOW OF CHARTS

FIGURE 73. FLOWVEAIT OF SUBROUTIHE 'GRAPHED'

 

 

 

 

FCR D FFEREKT
LOCATIONS:
T],T2,T3,T4,T§

1

‘DRAW STARTS OF I

TELP. EISTORIES

\START A HEW GRA-|

 

 

PEIC PAGE

 

 

 

SET UP THE WIN- fr
Dow OF CHART

I

DRAW CHART OF
TEMP. DIFFEREN—
CE VS. TIME

1

CLOSE OUTPUT
GRAPHIC FILE

‘ RETURN >

FIGURE 73. (CONT'D.)

 

 

 

 

 

 

<:¥EITER:)

DIVIDE THE SYS-
TEN INTO RECTA-
GULAR RFCICLS
PIID SI.fiLL G.‘ RIBS
AND ASSIGN TELP.

TO EACH NODE
DRAYJ BASIC CON-
FIGURPHTIOI OF

THE SYSTEM
(SUB. GapEETI)

 

 

 

 

 

 

 

 

 

   

 

 
 

OF EACH ISOTH-
RA}: I'I : 1 g o 0 III;

FOR EACH SIU LL
GRID OI? UPPER

BOUNDARY:
I = I,..ND(K)-1

FIGURE 74. FLO HotART OF 8"“DUTITE 'ISCTEERR'

 

 

Ox
C“

FIRD TRE START-
IRO POINT OF
ZSOTHERH ASE
323w THE ISOJH-
ERR TRROUOR THE
WHOLE REGION
(SUB. SIDET,SI-
DEB,SIDEL,SIDER)

 

 

 

 

 

 

 

 

FIRE THE START-
IRO POINT OF
ISOTHERH AND
DRAW THE ISOTH-
ERR THROUGH TH
WHOLE REGION
(SUB. SIDET,SI-
DEB,SIDEL,SIDER)

 

 

 

 

 

F OUR: 74. (COHT'D.)

 

  

:‘31{E;C:*I‘.SIL"1L
l\_‘I..I_J C’I: TJOLIBL‘IR

—-T\.\-v

_'\.I~.--LJAIA;\10

l : 1,..L.L}(I<)-1

 

1 _

 

FIND THE START-
ILC POI‘IT OF
ISO‘TI SEE-2&1 Alf-D
L‘RAL TEIE ISCTH-
EPE TEEROUG:' THE
'11: 7* " ’O“:

- ”'0".
I35. . L.) \UKJ- I\.4ILI

- ~:~ 0*, . "*‘r
('JVIIIJ. DlUZIiI,U-'~-

23,519RL ,SIDER)

 

 

 

 

 

s;R Racn SHALL
OR:O HO LRRR
Lwc ID:%1

J z 1,..ND(K)-1

 

l

 

 

FIND TrIE ST!‘ :RT-

'RAH TR: ISOTH—
ERH THROUGH THE
RROLL RROIOR

(SUB..SIDET,SI-
DEL,SIDEL,SIDER)

 

 

 

FIG’RE 74. COIL

T'D.

 

 

o

P.)

éda

Q\
07>

 

 

 

< 7 °.1r'*""fi
REL U11

TIGURE 7Q.

-
A!
.‘

(CONT'D)

 

 

 

 

169

ENTER

 
  
   
   
   
 
   
 
 

   

DOES
ISOTHERM PASSE
RIGHT BOUNDARY
OF THIS GRID?

   

  
 
   
 

NO
SET ITSELF YES DOES
EQUAL TO 1 ' ISOTRERM PASS
BECAUSE ' UPPER BOUNDARY
SUB. CAN OF THIS GRID°
NOT CALL
ITSELF

I NO

 

‘ RETURN ’

 
   

   
 
   
   
 

DOES
ISOTHERM PASSES
LEFT BOUNDARY
OF THIS GRID?

N O

 

QETURN ’

YES

YES

 

 

CONNECT ISOTH-
ERH SEGMENTS
AND CONTINUE
DRAWING ISOTH-
ERM FOR ADJAC-
ENT GRID

(SUB. SIDEL)

 

 

 

CONNECT ISOTH-

ERM SEGMENTS
AND CONTINUE
DRAWING ISOTH-
ERM FOR ADJAC-
ENT GRID '
(SUB. SIDER)

 

1

FIGURE .75. FLOWCHART OF SUBROUTINE 'SIDEB'

 

 

 

170

K-DIRECTIOI‘I OF
.1‘ BLOCK OF SILI-
ILER I‘IODES:
I = I1,..I2

DI RECTIOI‘I‘ CF

1351.. CCK OI“ SIM-

'1-
n
h
IL. 1R T’ “HA.
J

 

    
 
  

 

 

.'\\ U140 0

~ J1,..LT2
CPLLCJLiflL T: I I:C)R::-
KILLED TEX P.
UEJ(I,J) WIT}:

FORWARD-DIFFER-
ENCE EQUATION

 

 

 

 

 

 

 

 

‘ RETURN)

FIGURE ,76.FLOWCHART OF SUBROUTINE 'INNER'

 

 

 

 

IKKER, SBTK,....
(FORWARD-DIFFERENCE
EQUATION SUBROUTINES)

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

FIGURE 77. RELATIOZ'I AI-IOI‘IG PROGRAM UNITS

1 72
CDCCCCCCCC PRCBRAM DESCRIPTION CCCCC

THIS IS MAIN PRmRAM OF THE COMPUTER PRmRAM FOR 'COMPUTER
SIMILATION OF CRYOHICXOSCOPE CONDUCTION STAGE'

JCB :SERVE AS COMMANDING (ENTER OF THE WBCLE SET OF PRCORAM.

OOOOOOOOO
O00

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
CCCCCCCC VARIABLE IDENTIFICATION CCCCCCCC
II :LARGEST VALUE OF X-DIRRECTION SUBSCRIPT OF NODE NUDBR.
II :LARGEST VALUE OF Y-DIRRECI'ION SUBSCRIPT 0F NODE NUIBHI.
LX :DIIENSION IN X-DIRRECI'ION OF THE PLATE.
LY :DIIENSION IN Y-DIRRECTION OF THE PLATE.
DX :ARRAY TO STORE INTERVALS BETWEEN X-COORDINATES.
DY :ARRAY TO STORE INTERVALS BETWEEN Y-COORDINATES.
RT :ARRAY TO STORE TIME RATIOS.
RK :ARRAY TO STORE THERMAL CONDUCTIVITY RATIOS.
RR :ARRAY TO STORE THERMAL DIFFUSIVITY RATIOS.
BL :ARRAY T0 STORE BIOT NUIBERS FOR LEFT BGJNDARY.
BR :ARRAY T0 STORE BIOT NULBERS FOR RIGHT BWNDARY.

UO :2-D ARRAY TO STORE NORMALIZED TEMPERATURE DISTRIBUTION
OF FORMER TIME STEP(OLD ONE).

UN :2-D ARRAY TO STORE NORMALIZED TEMPERATURE DISTRIBUTION
OF PRESE‘JT' TIME STEP(NEW ONE).

RL :DIMENSION RATIO=LXI LY.
RTJIT :(REFER DUMMY ARGUMENTS OF SUBRGJT‘INE 'PLT. ISO.$')
BBC :BIOT NUIBER FOR BOTTOM SIDE OF COPPER.

BTCC :BIOT NUDBER FOR FORCED CONVECTION BY CO(LANT.

GOOOCOOOOOOOGOOOOOOOOOOOOOOGOOOOOOOOOOQO

BIO :BIOT NUDER FOR TOP SURFACE OF HEATING GLASS.

173

c
c BBQO :BIOT NBmBR FOR BOLB REGION OF Bomu SIDE OF
c GLASS BBATBR.
C
c BmAO :BIOT NUIBER FOR CONDUCTING BPBBcr BB'mBm SILVER
c EPOXY AND ITS CONTACTING COPPER.
C
C c
c c
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
csccccccc my _S'IORAGE BLOCK BLOCK 0000
c
PROGRAM N3PL'r
PARAMETER II-11.JJ=10.NM=10.M=5.IOI=10
REAL LX.LY.K(JJ—1).RK(JJ).RR(JJ) .RT(JJ).DX(II).R(JJ-1) .
_ DY(JJ) .BLUJ) .BRUJ)
BmaNSION 11(NN).12(MM).J1(NN).12(NN)
COMMON /C1/TNAK.TN.UN(II.JJ)
_ /C2/0T.RL.00(II.JJ)
_ lC3/VX.VY,VZ.RX.RY,RZ
_ IC4/K.R
_ IC9/0,BmIT.0'r.uB,0c.B'r,m,m..BR.BC,BTAO.m0
_ /c5/rom=, no, IDA’I‘A, IOFF, 13m. ICDN'IOUR, IZDG. INIT
__ lC6/L.G'I'EM(8.250)
_ lC7/XO.X1,Y0.11,NX,NY.XZ.XS,Y2,Y3,IX.HY
_ Ice/Gnumzsm
_ lC10/X(II).Y(JJ)
_ /C21/XC(6.NK).YC(6.NK).ZC(6.6.NK)
_ IC22/BPBI.PBI(NN).NB(NK)
_ /c23/ ITSBLF. IT. JT.KT. nn'
_ lC24/TEM(II.JJ) .11:an
_ lC25/LX.LY
_ lC26/II3DRG.11,IZ,11.JZ
c c
c C
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
C c
CICCCCCCCC INITIALIZATION BLOCK BLOCK 0100
c c
C c
CCCcccccccccccccccccccccccccCCccccccCCccccccccccccccccccccccccccccc
c C
CPCCCCCCCC PROCESS BLOCK BLOCK 0200
C
C----DA'I‘A INPUT
c
CALL INOB'r
c

C-----INITIALIZING THE GRAPHIC MOE

C
c-
C

20

C
c-
C

14

19

16

17

 

 

171+

IF (I3DG.NE.O) THEN
CALL GRAPBBD
ELSE
CALL GRAPBZD
END IF

IN ITIALIZATION

TM=O.
I(OUNT=O
IF].
IF (INIT.m.O) THEN
DO 1 I=1.II
D0 1 181,11
UO(I.J)=UINIT
CONTINUE
ELSE
DO 20 I=1.II
DO 20 J’=1.IJ'
UO(I.J)=(TEN(I.J)+200.)[250.
CONTINUE
BUD IF
TEIDLD=GTEM(1.1)
CALL (RITFTEM

CONVERTING INPUT PARAMETERS INTO REQUIRED FORMS

RIFLX/LY
DO 14 J=1.JJ'-1
RT(J)=T|D‘R(J)/LX/LX
DY(J)=Y(J+1)-Y(J)
CONTINUE
DO 19 I=1,II-1
DX(I)=X(I+1)-X(I)
CONTINUE
DO 16 182.11-1
RK(J)=K(J') /K(J-1)
RR(J)-R(J)lR(J-1)
CONTINUE
BBC==BB‘LX/K(1)
BTCC=BC‘LX/K (1)
BTO=BT‘LX/K ( 5)
Bm0=-O‘LX/K ( 5)
BTOAG=BTAG’LXIK(5)
DO 17 186.10
BL( J) =BL‘LX/K (J—l)
CONTINUE
DO 18 182.5
BRU) =RR‘LX/K (J-l)

175

18 CONTINUE

 

C
C----CALCULATION OF TEMPERATURE DISTRIBUTION WITH FORWARD
C- DIFFERENCE METHOD
C
13 'TMETM+DT
ICOUNT=ICOUNT41
C
C----CONSTRUCTING FORWARD DIFFERENCE EQUATIONS
C
I81
DO 2 182.6

CALL SBTHCRT(J).DX(I-1).DX(I).DY(J).I.I.J.1..1..0..UB.BBC)
2 CONTINUE
J=5
DO 3 188,10
CALL SBTHLRT(J).DX(I-1).DX(I).DY(J).I.I.J.1..1..0..UB.
__ BBQO)
CONTINUE
1810
DO 4 186.10
CALL SUP(RT(J-1) DDX(I-1) DDX(I)JDY(J.1) '1' 1.1.1. :1. no. .UT,
_ BID)
4 CONTINUE
CALL SUP(RT(1).DX(1).DX(2),DY(1),2.2.2.1..1..0..UC,
fiTCC‘DX(1)/(DX(1)+DX(2)))
CALL SUP(RT(5).DX(3).DX(4).DY(5).4.4.6.1..1..0..UT.BTQAG)
I=3
DO 5 J=3.5
CALL SLFr(RT(J-1).DX(I).DY(J-1).DYU).I.J.J.RK(J).RR(J).
_ O..O..O.)
5 CONTTNUE
I=5
DO 6 187.9
CALL SLFT(RTTJ-1).DX(I).DY(J-1).DY(J).I.J.J.RK(J).RR(J).
.. 0..UT.BL(J))
6 CONTTNUE
1811
DO 7 1:7,9
CALL SRGT(RT(I-1).DX(I-1).DY(J-1).DY(J).I.J.J.RK(J).RR(J).
_ 0..O.,O.)
7 CONTTNUE
187
D0 8 J=2,4
CALL SRGT(RT(I-1).DX(I-1).DY(J-1).DY(J).I.J.J.RK(J).RR(J).
_ 0..IIB.BR(J))
8 CONTTNUE
CALL SRGTTRT(5).DX(10).DY(5).DY(6).11.6.6.RK(6).RR(6).
§!DX(10)/2..0..0.)
DO 9 1:2.5

w

1'76

DO 9 I=4.6
CALL INNER(RT(J-1).DX(I-1).DX(I).DY(J—1).DY(J). I. I.J. J.
_ RK(J).RR(J>.0.)
9 CONTINUE
no 10 J=7.9
00 10 I=6.10
CALL INNER(RT(J-1).DX(I-1) .DX(I).DY(J-1) .DYU) . I. I. J. J.
_ RK(J).RR(J).0.)
10 CONTINUE
J=6
00 11 I=6.10
CALL INNER(RT(J-1).nx(I-1).nx(1).BT(J-1) .DYU). I. I.J.J.
_ RKU).RRU).G‘(DX(I~1)+DX(I))12.)
11 CONTINUE
CALL CLU(RT(1).DX(1).DY(1) .1.2.0. .UT.UC.BR(2) .BTCC)
CALL CLU(RT(5).DX(3).DY(5).3.6.0..0..UT.0..B‘IQAG)
CALL (1.0(RT(9) .DX(5).DY(9) .5.10.0..UT.UT.BL(10) .8110)
CALL CRD(RT(1) .DX(6) .DY(1).7.1.0..IB.IB.BR(2) .BBC)
CALL CRD(RT(5) .DX(10) .DY(5) .11.5.0. .0..1m.0. .BBOO)
CALL (LNRTU).DX(1).DY(1).1.1.0..UT.UB.BR(2).BBC)
CALL CRU<RT(9).0x(10).mmn.11.10.0..0..UT.0..Bm)
CALL ICLU(RT(5).DX(4).DX(5).DY(5).DY(6).5.6.RK(6).RR(6).
§‘(DX(4)+DX(5) )IZ . .UT.UT.BL( 6) mm)
CALL ICLU<RT<1>.0x<2).0x<3).BT(1).0Y(2).3.2.RK(2).RB(2).

 

Q..O..O..O..O.)
CALL ICRD(RT(4) .DX(6) .DX(7) .DY(4) .DY(5) .7 .5 .RK(5) .RR(5) .0. .
EB.IB.BR(5).BBQO)
C
C- STORE CALCULATED TEMPERATURE DISTRIBUTION INTO 'UO'
C
DO 15 I=1.II
DO 15 J=1.JJ
UO(I.J)=UN(I.J)
15 CONTINUE
C
C---—CBECK WHETHER IT IS THE THE TO CUT BEAT GENERATION OFF
C-----OR WTPUT THE DATA
C
IF (TN.GT.TOFF) THE]
G=O.
IF ((IWUNI'IIOFF)‘IOFF.NE.I(OUNT) GO TO 13
ELSE
IF ((ICOUNT/IDATA)‘IDATA.NE. MOUNT) GO TO 13
BID IF
C
C----WTPUT TEE TEMPERATURE DIS'mIBUTION
C
IFL-Pl
CALL (OTPI‘EM

177

C----OUTPUT THE GRAPES EITHER THREE DIMENSIONAL TEMPERATURE
C-----PROFILE OR ISOTHERM ON “0 DILENSION
C

IF (I3DG.NE.O) THEN
IF ((ICOUNT/ISDG)‘I3DG.m.I(OUNT) THEN
CALL NEWPAG
CALL GSDENTT
CALL HOME
CALL ANBDDE
END IF
ELSE IF (I(ONTOUR.NE.O) THEN
IF ((HDUNT/IGWTOURPICONTOULFD. ICDUNT) THEN

CALL NEWPAG
CALL ISOTHERM
CALL HOME
CALL ANIDDE
END IF
END IF
C
C------CHECK WHETHER THE TILE IS UP
C
IF (TM.LT.TMAX) GO TO 13
C
C----WTPUT 'EE TEMPERATURE VS. TIME PLOTT‘DNG AND TEM DIFFERENCE
C-----VS. TIME PLOTT'ING
C
LL=L
IF (I3DG.NE.0) CALL GRAPH2D
IF (12DG.FD.1) THEN
CALL NEWPAG
CALL GZDENTTUL)
END IF
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CECCCCCCCC ERROR MESSAGES BLOCK 0900
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CFCCCCCCCC FORMAT STATEMENTS BLOCK 9000
C C
C C

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C

CALL EXIT

END

178

CDCCCCC SUBRGJTINE DESCRIPTION CCCCC
C
C JIB :TAKE CARE OF DATA INPUT AND OUTPUT.
C C
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
CCCCCCCC VARIABLE IDENTIFICATION CCCCCCCC

VX.VY.VZ :COORDINATES OF VIEWPOINT FOR 3-D GRAIH.
RX. RY.” :COORDINATES OF REFERENCE POINT FOR 3-D GRAPH.
M3DRG :NUIBER OF REiION SYSTEM DIVIDED FOR 3-D PLOT'TING.

11.12 :ARRAY STORE X-DIRRPL'TION SUBSCRIPTS OF TWO LIMITING
NODES OF A PARTICULAR REIGN.

11,12 :ARRAY STORE Y-DIRRECTION SUBSCRIPTS OF TWO LIMITING
NODES OF A PARTICULAR REEION.

K :ARRAY STORE NUMERICAL VALUE OF CONDUCTIVITY OF EACH LAYER.
SUCH AS K(1) IS THAT OF COPPER.

R :ARRAY STORE NUMERICAL VALUE OF DIFFUSIVITY OF EACH LAYER.
SUCH AS RU) IS THAT OF COPPER.

G :RATE OF HEAT GENERATION.

UINIT :INITIAL TEMPERATURE OF THE SYSTEM.(IF HODDGENEOUS)
UT :ADBIENT TEMPERATURE ON TOP OF THE SYSTEM.

UB :ADBIENT TEMPERATURE BENEATH THE SYSTEM.

UC :COOLANT TEMPERATURE.

HT :HEAT TRANSFER COEFFICIENT FOR TOP SURFACE OF THE SYSTEM.
HB :HEAT TRANSFER COEF. FOR BOTTOM SURFACE OF THE SYSTEM.
HL :HEAT TRANSFER COEF. FOR LEFT SIDE OF THE SYSTEM.

HR :HEAT TRANSFER COEF. FOR RIGHT SIDE OF THE SYSTEM.

BC :HEAT TRANSFER COEF. FOR FORCED CONVECTION BY COQANT.

HTAG :HEAT TANSFER COEF. ACCOUNT FOR THE (ONDUCTING EFFECT
BETWEm SILVER EPOXY AND ITS CONTACTING (OPPER.

GOGOOOOOOOOOOOOOOOOOOOOOCOOCOOOOOOOGOGOOGQO

OOOOOOOOOOGOOOOOOGOGOOGOGOOOOOGOOOOOOG000060600600

179

HBO :HEAT TRANSFER COEF. FOR THE HOLE ON THE BOTTOM SURFACE.

DT :TIME DIFFERE‘JCE FOR EACH TIME STEP ADVANCED IN
CALCULATION WITH FORWARD DIFFERENCE METHOD.

TOFF :TIME TO CUT HEAT GWERATION OFF.

TMAX :TOTAL TIDE DURING WHICH THE RESPONSE OF THE SYSTEM
WANT TO BE PREDICTED.

TDD :ARBITRARY VALUE BASE ON WHICH TIME IS NORMALIZED.

IDATA :CONTRmLING INTHEER FOR DATA OUTPUT. IT WILL (IITPUT
THE DATA EVERY 'IDATA' TIIB STEPS OF CALCULATIONS.

IOFF :CONTRmLING INTHEER FOR DATA OUTPUT AFTER HEAT BEING
CUT OFF. IT WILL (DTPUT HE DATA EVERY 'IOFF' TIDE
STEPS OF CALCULATIONS.

I3DG :CONTRmLING INTHSER FOR THREE DIMENSIONAL GRAPH OUTPUT
IT WILL HOT 3-D GRAPH EVERY 'ISDG' TIME STEPS OF
CALCULATION. IF IBDGBO. STOP THE (RJTPUT.

ImNTOUR :CONTRmLING INTEGER FOR CONTOUR LINEUSOTHERM)
OUTPUT. IT WILL PLOT ISOTHERM EVERY 'ICDNTOUR' TIME
STEPS. IF ICONTOUEO. STOP THE OUTPUT.
(3-D GRAPH AND ISOTHERM CANNOT BE REQUIRED AT THE SAME TIME)
I2DG :CONTRmLING INTHiER FOR CHART OUTPUTCTEM. VS. TIDE AND
TEM. DIFF. VS. TIME). IF 12DG=1, DEMAND THE OUTPUT.
IF 12DG=0. STOP THE OUTPUT.
INIT :CONTRmLING INTEGER TO INDICATE WHETHER THE INITIAL
TEMPERATURE DISTRIBUTION IS HODDGENIEOUS OR NOT.
IF IT IS HODDGEVEWS, THE] INITBO. IF NOT, THE INIT‘=1.
X :ARRAY STORE X-COORDINATES.
Y :ARRAY STORE Y-COORDINATES.

XO.X1.YO,Y1 :MINIDIJM AND MAXIIIIM VALUES FOR FRAME OF
TEMPERATURE VS. TIME CHART.

X2.X3.Y2.Y3 :MINIMUM AND MAXIIIJM VALUES FOR FRAME OF
TEM. DIFFEREWB VS. TIME CHART.

NLNY :NUDBER OF GRIDS DIVIDED BETWEEN TWO LIMITS OF TEM.
VS. TIDE CHART.

180

MX.MY :NUEER OF GRIDS DIVIDED RENEW TWO LIMITS OF TEM.
DIFF. VS. TIE CHART.

TEM :2-D ARRAY TO STORE TEMIERATURE DISTRIBUTION OF SYSTEM.
(IN DEGREE CELSIUS)

GTEM. GTIE :ARRAY TO STORE TEMPERATURE AND CORRESPONDING TIE
FOR LATER PLOTTINGS OF TEM. VS. TIE CHART AND
TEM. DIFF. VS. TIE CHART.

L :INTEGER FOR GTEM AND GTIE TO INDICATE WHICH TIE STEP IS
GOING RIGHT NW.

TEDDLD :STORE TEMPERATURE OF FORMER TIME STEP FOR RATE OF
TEM. CHANGE CALCULATION.

OOOGOOOOOOOOOOOOOOOOOOG

c
c
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c
sccccccc mm _smma BLOCK BLOCK 0000
c
sunnourmn moor
pnmm II=11.JJ=10.NM=10.NK=5.MM=10
REAL K(JJ-1).R(JJ-1).LX.LY
DIMENSION n(uu).n(uu).11(uu).12(um
comm lC1/TMAX.TM.UN(II.JJ)
_ IC2/DT.RL.UO(II.JJ)
_ lea/vx.vr.vz,nx,u,nz
_ lC4/K.R
_ I C5/TOFF. Tm. IDATA. 10131:. 1300, 1001010011. 1200. INIT
_ lC6/L,GI'EM(8.250)
__ IC7/XO,X1,YO,Y1,NX,NY,XZ.X3.Y2,Y3.MX.MY
_ Ics/G'rnmzsm
_ lC9/G,UINIT,UT.UB.UC,HT,B,HL,HR, nc,1ms,m0
_ lClO/X(II).Y(JJ)
_ /C22/DPHI,PBI(NM),ND(NK)
_ lCZS/LX.LY
_ lCZ4/TEM(II.JJ) .113.an
__ Iczsmanne. n .12 .11,12
0 c
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c c
cmccccccc INITIALIZATION BLOCK BLOCK 0100
c c
c c
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c c

CPCCCCCCCC PROCESS BLOCK BLOCK 0200

181

OPEN(11.FILE='IM3PLT")
OPEN(10,FILE='OM3PLT')

C

C----—READING THE INPUT DATA

C
READ(11,‘)VX.VY.VZ.RX.RY,RZ
READ(11.‘)M3DRG.(11(M).IZ(M).Jl(M).12(M).M=1.M30RG)
READ(11.‘)K.R
READ(11,‘)G,UINIT.UT,IB,UC,HT,B.HL.HR.HC.HTAG.BO
READ(11,‘)DT.TOFF,TMAX.TMO
READ(11.‘)IDATA.IOFF.I3DG.ICONTOUR.IzDG.INIT
READ(11.‘)(X(I).I=1.II)
READ(11.‘)(Y(J).J-1.JJ)
READ(11,‘)XO,X1,Y0.Y1,NX,NY.XZ.13.Y2.Y3,MX.MY
READ(11.‘)(PHI(I).I=1.NM)
IF (INIT.NE.0) THEJ

DO 4 J=JJ.1,-1
READ(11.‘)(TEM(I.J),I=1,II)

4 CONTINUE
END IF

C

C----OUTPUT THE DATA REID FOR KEEPING RECORD ON COMPUTER OUTPUT

C
WRITE(10.'(/. " VIEWPOINT AND REFERENCEPOINT ARE' './)')
WRITE(10,‘)VX,VY,VZ,RX,RY,RZ
WRITE(10,'(/."M3DRG II 12 11 12 ARE",/)')
WRITE(10,‘)M3DRG.(11(M),12(M).JI(M).J2(M),M#1.M3DRG)
WRITE(10.'(/."K(JJ) AND R(JJ) ARE",/)')
WRITE(10.‘)K.R
WRITE(10.'(/." G UINIT UT'UB UC HT HB HL HR.HC HTAG HBO ARE

_ ' './) ')

WRITEUO,’)G,UINIT,UT,UB,UC,HT,HB,HL,HR.HC,HTAG.EO
WRITE(10.'(/." DT TOFF TMAX TDD ARE"./)')
WRITE(10,’)DT,TOFF,TMAX,TMO
WRITE(10.'(/." IDATA,IOFF.13DG.ICONTOUR,IZDG,INIT ARE"./)')
WRITE(10,‘)IDATA.IOFF.IBDG,ICONTOUR.I2DG.INIT
WRITE(10.'(/."X COORDINATES ARE"./)')
WRITE(10.‘)X
WRITE(10.'(/."Y COORDINATES ARE"./)')
WRITE(10.’)Y
WRITE(10.'(/."XO X1 Y0 Y1 NX.NY X2 X3 Y2 Y3 MX MY ARE"./)')
WRITE(10,‘)XO.XI,YO,Y1,NX.NY,XZ,X3,Y2.Y3,MX,MY
WRITE(10.'(/."PHI(NM) ARE".I)')
WRITE(10.‘)(PHI(I).I=1.NM)

C

 

NORMALIZING THE,COORDINATE'VALUES

LX=X(II)
LY=Y(JJ)

 

DO 2 I=1,II
X(I)=X(I)ILX
2 CONTINUE
DO 3 181,11
Y(J)=Y(J)/LY
3 CONTINUE
CLOSE(11)
RETURN
C
C SECOND ENTRY POINT FOR TEMPERATURE DISTRIBUTION WTPUT
C
ENTRY OUTPTEM
C

C----CHANGE NORMALIZED TEMPERATURE TO CELSIUS DEGREE
C
DO 21 181.11
DO 21 J=1.JJ
TEM(I.J)=250.‘UO(I,J’)-200.
21 CONTINUE

 

RTIFTM'TK)
C
C-----STORE TEMPERATURES OF CERTAIN LOCATIONS FOR LATER PLOT'TING
C
GTEM(1.L)=TEM(II.J'J)
GTEM(2.L)=TEM(II-4.J’J)
GTEM(3.L)=TEM(II-6.JJ)
GTEM(4.L)==TEM(7.1)
GTEM(5.L)=TEM(11.5)
GTIE(L)=RTM
C
C OUTPUT THE TEMPERATURE DISTRIBUTION
C
WRITE(1.‘)RTM,TEM(II,IJ') .TEM(II-4,JJ) .TEM(II-6,J'J) .TEM(7.1) .
__ TEM(11.5)
WRITE(10.202)RTM

202 FORMATU. 'THE TEMPERATURE DISTRIBUTION IS AT TIME-“.FBA)
DO 1 J-JJ',1.-1
WRITE(10.201) (TEM(I.J) . I=1. II)

201 FORMAT(3X.11(F7.1))

1 CONTINUE

C

C------OUTPUT THE RATE OF TEMIERATURE CHANGE

C
TEEATE=(TEM( II. II) -TEDDLD) ITID/DT/IDATA
TEIDLD=TEM(II.JJ)
WRITE(10.203)TEEATE

203 FORMATU. 'RATE OF TEMPERATURE CHANGE IS' .2X.F14.2.1X. 'C/SHZ')
IF (TM.LT.TMAX) RETURN
CLOSEUO)

183

C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CERCCCCCCC ERROR MESSAGES BLOCK 0900
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CFCCCCCCCC FORMAT STATEMENTS BLOCK 9000
C C
C C

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C

RETURN

END

184

 

CDCCCCC SDBRODTINE DESCRIPTION CCCCC
c
C ma :(1)ACI'UATE THREE DIMENSIONAL GRAPHIC MODE.
C
c (2)PL0T 3-D TEMPERATURE DISTRIBUTION PROFILE.
C
C C
C C
c c
CCCCCCCCCCCCCCCCCCccccccccccccccccccccCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CVCCCCCCCC VARIABLE IDENTIFICATION CCCCCCCC
C
C m :MAXIMUM NUIBER OF REGIONS SYSTEM CAN BE DIVIDED FOR
C 3-D PLO’l'I‘mG.(CAN BE CHOOSEN FREE!)
C
c c
C C
CCCCCCCcccCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCcccccccc
c C
CSCCCCCCC ENTRY _STORAGE BLOCK BLOCK 0000
C C
C
C------FIRST my POINT GRAPH3D Is USE) TO ACTUATE GRAPHIC MODE
C mm
c
SUBRODTINE GRAPBSD
PARAMETER II-ll.JJ=10.MM=10
CBARACTER‘3O LABTM
REAL K(JI-1).R(JJ—1)
DIMENSION I1(NN).Iz(MN).n(uM).Jz(uM)
COMMON lC1/TMAX.TM.UN(II.JJ)
_ /C3/VX,VY.VZ.RX,RY,1E
_ IC4/K.R
_ IC6/L.GTEM(8.250)
_ lC8/GTIE(250)
_ 109/0.nmIT.uT,m.DC,BT,m,RL.BR,BC,BTAG.BBO
_ IC10/X(II).Y(JJ)
_ IC26/M3DRG. 11 . 12 .11 .12
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
c c
CICCCCCCCC INITIALIZATION BLOCK BLOCK 0100
C C
c c
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CPCCCCCCCC PROCESS BLOCK BLOCK 0200
C

CALL INITT3 ( 120)

185

CALL 0PEN’I'K('63D'.I)
CALL CAR'IVP(VX.VY.VZ)
CALL LOOKAT(RX.RY.RZ)
CALL HAGNFY(16.)

RETURN
C
C-----SECOND ENTRY POINT G3DENT1 IS USE TO PLOT 3-D TEMPERATURE
C----PROFILE
C
EJTRY GBDEN‘I‘I
C
C—----DRAW INC THE REFEREE FRADE
C

CALL IDVEAB(O..0..0.)
CALL VECTOR(1.2.0..O.)
CALL mVEA3(0..0..0.)
CALL VECTOR“). .1.2.0.)
CALL IDVEA3(0..0..0.)
CALL VECTOR(0..0. .13)
CALL KWEA3(1.3.0..0.)
CALL (BAR'IX( '1' .0.7)
CALL mVEA3(O..1.3.0.)
CALL GARH('Y' .0.7)
CALL IDVEA3(O..0..1.3)
CALL CHAR'IX('TEM'.0.7)
CALL mVEA3(0.013 .0..0.5)
CALL DRAWA3(-0.013 .0. .0 .5)
CALL mVEA3(-0.1.0..0.5)
CALL CHARTK( '0.5'.0.7)
CALL mVEA3(0.013.0..1.)
CALL DRAWA3(’0.013.0..1.)
CALL mVEA3(-0.1.0..1.)
CALL CBARTK('1.0'.0.7)
CALL IDVEA3(0.5 .0.."0.013)
CALL DRAWA3(0.5 .0..0.013)
CALL IDVEA3(0.5.0..-0.1)
CALL GARH('0.5'.0.7)
CALL mVEA3(1..0..0.013)
CALL DRAW” (1 o .0 O .‘0 0013)
CALL mVEA3(1..0..’0.1)
CALL (HARK('1.0'.0.7)
CALL IDVEA3(0.013.0.5 .0.)
CALL DRAWA3(-0.013.0.5 .0.)
CALL mVEA3(’0.1.0.5.0.)
CALL (BAR'I'K('0.5'.0.7)
CALL DDVEA3(0.013.1..0.)
CALL DRAWA3(-0.013.1..0.)
CALL mVEA3(-0.1.1..0.)
CALL (EARTK('1.0'.0.7)
WRITE<LABT3L101J GTIIE(L)

 

on/
ICC

101 FORMAT(2X,'TIIE IS '.F6.1,' SEC')
CALL mVEA3(1.1.0..O.3)
CALL CHARTKUABTILOJ)
C
C-“u-PLOTI'ING THE TEMPERATURE DIS'mIBUTION HOFILE
C
DO 7 N=1.NBDRG
DO 3 J-Jl(M).J'2(l)
CALL IDVEA3(X(II(N)).Y(J).UN(II(I).J))
DO 3 I=II(M).IZ(M)
CALL DRAWA3(X(I).Y(J).UN(I.JH
3 CONTINUE
DO 4 I'II(N).IZ(')
CALL DVEAMXU).Y(Jl(l)).UN(I.Jl(l)H
CALL DRAWA3(X(I).Y(J).UN(I.J))

4 CONTINUE
7 CONTINUE
CALL AN WDE

IF (TN.LT.TIAX) RETURN
CALL C'LOSTK(I)

C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CERCCCCCCC ERROR MESSAGES BIDCK 0900
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CFCCCCCCCC FORMAT STATEMENTS BLOCK 9000
C C
C C

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C

RETURN

END

187

CDCCCCC SUBROUTINE DESCRIPTION CCCCC
c
C ms :(1)ACTUATE 'NO DITENSIONAL GRAPHIC MODE.
C
C (2)PLOT ISOTHERM.
C
C (3)PLOT TEMPERATURE vs. TIME CHART.
C
C (4)PL0T TEMPERATURE DIFFEREKE vs. TIME CHART.
C
C c
c C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCccccCCCCccccCCCCCCCCCCCCCCCCCCCCCCCCCCCC
c C
CVCCCCCCCC VARIABLE IDENTIFICATION CCCCCCCC
C
C LABEL.LABTII :ARRAY 0F CHARACTER VARIABLES STORE NUMERICAL
C VALUES SO THAT IT CAN BE DRAWN 0N GRAPH.
C
C TEMLCI' :CHARACTER VARIABLE STORE NUMERICAL VALUE STAND FOR
C TEMPERATURE LOCATION .
C
C RGK :RANGE OF x-COORDINATE TO BE EXPRESSED m CBAR'I‘BXl—XO.
C
C RGY :RANGE OF Y-COORDINATE TO BE EXPRESSED IN CHARTaTI-TO.
C
C P.O :ONLY Two LOCAL VARIABLES USED IN D0 LOOP.
C
C C
ccCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
c C
CSCCCCCCC ENTRY _STORAGE BLOCK BLOCK 0000
C C
c
C----—-SUBROUTINE GRAPH2D IS US. TO ACTUATE GRAIHIC MODE ONLY
C
SUBROUTINE GRAPH2D
PARAMETER II=11.JJ=10
CHARACTER-10 LABEL(4.21) .IABTM.TEMLCI'60
COMMON /C6/L.GTEM(8.250)
_ /CT/K0.XI.YO.II.NT,NI.K2.13.I2.Is.Mx.MT
_ lCS/GTIME(250)
/C9/G,UINIT.UT.UB.UC.HT.IE.HL.HR.HC.HTAG.IBO
- lClO/X(II).Y(JJ)
c C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
c C
CICCCCCCCC INITIALIZATION BLOCK BLOCK 0100

C C

C

188

C

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC

C

CPCCCCCCCC

C

C
C

 

CALL
CALL

PROCESS BLOCK

INIT'I'(120)
0PENTK( '020' . I)

RETURN

C----CROSS-SECI‘ ION OF COG. ING STAGE

C

102

C
C

C

C

 

ENTRY GZDENTT ( LL)

CALL
CALL
CALL
CALL
CALL
CALL
CALL
CALL
CALL
CALL
CALL
CALL
CALL
CALL
CALL
CALL
CALL
CALL

D'INDO(-0 .75.2 35.-3.75.5 .25)
IDVEA(X(1).Y(1))
DRAWMXU) .Y(1))
DRAWA(X(7) .Y(4))
DRAWA(X(3).Y(4))
DRAWA(X(3) .Y(2))
DRAWA(X(1).Y(2))
DRAWA(X(1) .Y(1))
DDVEA(X(3) .Y(5))
DRAWA(X(11) .Y(5))
DRAWA(X(11) .Y(10))
DRAWA(X(5) .Y(10))
DRAWA(X(5).Y(6))
DRAWA(X(3) .Y(6))
DRAWA(X(3).Y(S))
mVEA(0.375.-1.0)
CEARTKUTIIE IS'.0.7)
IDVEA(0.625 .‘1.0)

WRITE(LABTM, 102) GTDBUJ.)
FORMAT(F5 .1)

CALL
CALL
CALL

CHARTK(LABTM.0.7)
IDVEA(0.85 .-1 .0)
(EARTK( 'SEC' .0 .7)

RETURN

ENTRY GZDENT2 (LL)

C
BLOCK 0200

SECOND ENTRY POINT GZDENTI FOR PLOTTING TEE wNTOUR OF

THIRD ENTRY POINT GZDENTZ IS USE) FOR PLOT'I‘INT TEMPERATURE
C----VS. TILE CHARTS AND TEM DIFFERENCE VS. TIME CHART

C-----PLOTTING TEMPERATURE VS. TIDE CHARTS FROM T1 T0 T5

C

RGX=X1-XO

RGY=Y1-Y0
DO 11 J=1 .5
CALL WINDOW .5‘ (Kl-+10) -MX. 0 .5 ‘ (ll-+10) +mx.0 .5‘ (Y1+YO)
-RGY.0.5‘(Y1+YO)+RGY)

189

C
C----DRAW BORIZCNTAL AND VERTICAL GRID LINES
C
DO 3 I=1 . NX+1
P=(I-1) ‘RGX/NX‘KO
CALL DOVEAU’. Y0)
CALL DRAWA(P. I1)
3 CONTINUE
DO 4 I=1.NY+1
Q=(I-1) ‘RGY/NY-t-YO
CALL IDVEAUO . Q)
CALL DRAWAUI . Q)

4 CONTINUE

C

C----WRIT ING CRARAT'ERS
C

CALL IDVEA( (11+XO) I2 .-mX/6. .Y1+mY/9 .)
WRITE(TEM.CI‘.103) J
103 FORMAT( 'TEMPERATURE ' .12 . ' VS TIME')
CALL (BARTK (TEMLCI'. 0 .7)
CALL DDVEM (X1+XO) l2 .-mX‘0 .875 .(Y1+Y0) I2 .)
CALL CRARTK( '(DEGREE C) ' .0.7)
CALL 101/EM (11+XO)/2 .-mX/24. .YO-O .222‘RGY)
CALL CEARTK( ' ( SEC) ' .0.7)
DO 7 I=1.NX+1.2
P=(I-1 . ) ‘RGX/NX+X0-mX/15 .
Q=(I-1 .)‘RGX/NX+X0
CALL IOVEA(P.Y0-mY/9.)
WRITE(LABH.( 1 . I) .101) 0
CALL (BARTKUABEJI . I) .0 .7 )
7 CONTINUE
DO 8 I=1.NY+1.2
PHI-1) ‘RGY/NY-I-YO
CALL IDVEMXO-mX/G . .P)
WRITE(LABE..( 2 . I) .101) P
CALL CRARTK(LABE.(2 . I) .0 .7)

8 CONTINUE
C
C----PLOTTING TEE CURVES
C

CALL IDVEA(GTIME(1) .GTEM(J.1))

DO 1 I=2.LL

CALL DRAWA(GTI|E(I).GT'EM(J.I))

1 CONTINUE

CALL NEWPAG
11 CONTINUE
C
C-----PLOTTING TEM DIFFERWCE VS. TIME CHART
C

CALL DWINDO(0.5‘(13+X2)-(13-12) .0 .5‘(13+XZ)+(13-X2) .

190

0.5. (Y3+Y2)-(Y3-Y2) .0 .5‘ (Y3+YZ) +(Y3-Y2))
C
C------DRAW HORIZG‘ITAL AND VERTICAL GRID LINES
C
D0 5 I=1.MX+1
P=(I-1) ‘ (KB-12) IMX+X2
CALL DVEAU’. Y2)
CALL DRAWA(P.Y3)
5 (DNTINUE
DO 6 I=1.MY+1
Q=( I-1) ‘ (Y3-Y2) IMY+Y2
CALL IDVEAUZ . Q)
CALL DRAWA(X3 . Q)

6 CONTINUE

C

C----WRITING TEE CHARAT'ERS
C

CALL DVEA( (X3+X2) I2 .-(X3-X2) I6. .Y3+(Y3-Y2) I9.)
CALL (EARTKUTEM DIFFERENCE VS TIME'.0.7)
CALL mVEM (KB-+12) I2 .-(X3-X2) ‘0 .875 . (Y3+Y2) I2 .)
CALL CHARTK('(DEGREE C) '.0.7)
CALL DVEM (13+X2) I2 .-(X3-X2) I24 . .Y2-0 .222‘ (YB-Y2”
CALL CRARTK( ’ ( SEC) ' .0 .7)
DO 9 I=1.MX+1.2
ku-l .)‘(X3-12) /MX+X2-(X3-X2) I15 .
Q=( I-1 . ) ’ (XS-12) IMX+X2
CALL DVEA(P.Y2-(Y3-Y2)I9.)
WRITE(LABH.(3.I).101) Q
CALL GART'K(LABEL(3.I).0.7)
9 CONTINUE
D0 10 I=1.MY+1,2
P=(I-1) ‘ (Y3-Y2) IIIY+Y2
CALL IDVEA(X2-(X3-X2) I6. .P)
WRITE(LABH.(4.I) .101) P
CALL CHARTX(LABE.(4.I).O.7)
10 CONTINUE
101 FORMAT(F6 .0)

C
C----PLOT'I'ING CURVES
C
CALL IDVEAULO.)
D0 2 I=1,LL
DTEM=GTEM(1.I)-GTEM(2.I)
CALL DRAWA(GTIIE(I) .UI‘EM)
2 CONTINUE
CALL C'LOSTKU)
CALL ANIDDE
C C
C C

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC

191

C C
CERCCCCCCC ERROR MESSAGES BLOCK 0900
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CFCCCCCCCC FORMAT STATEMENTS BLOCK 9000
C C
C C

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C

RETURN

END

1 92
CDCCCCC SUBRWTINE DESCRIPTION CCCCC

JCB :DRAW ISOTEERM OF DIFFERET T'EMERATURE ON A PRESCRIBED
REIGN.

00006600
600

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
CCCCCCCC VARIABLE IDENTIFICATION CCCCCCCC

20

NM :TOTAL NUBER 0F REIGN WHICH WANT TO BE DRAWN ISGTBERM
0N THEM

NK :TOTAL NULBER OF ISOTEERM 0F DIFFERET TEMPERATURE

LABEL :CBARACI‘ER VARIABLE FOR LABflJNG TEMPERATURE OF
ISGTEERM.

XC,YC,ZC :COORDINATE AND FUNCTION VALUES 0F SELECTED NWES
FOR DRAWING ISGTBERM

DPHI :POSITIVE SMALL NUDBER FOR DECIDING WHETHER A STARTING
POINT IS BETWEE TWO NODES

PEI :TEMPERATURE VALUES 0F ISOTHERM
ND :GRID NUDBER OF EACH REIGN

ITSEF :CRITERIA VARIABLE FOR SUBRWTINE T0 CALLING ITSEF
BECAUSE SUBRGJTINE IN FORTRAN CAN'T CALL IT'SEF

IT.JT :REFERENICE NGDE NUDBER FOR CALLING A SUBRUTINE T0 DRAW
ISGTEERM:

FOR 'SIDEB'. THEY ARE DEFINED AS THE NODE NUDBER GF
LEFT-DOWNWARD CORNER OF A CERTAIN REIGN.

IN THE SALE WAY. 'SIDET' IS L-U ERNER.

'SIDEL' IS L-D CORNER.

'SIDER' IS R-D CORNER.
WHERE ' SIDEB ' , ' SIDET' , ' SIDEL' , ' SIDER' ARE SUBRCIJTINES
T0 DRAW ISOTIIERM.

C

C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC

C
SCCCCCCC ETRY _ST’GRAGE BLOCK BLOCK 0000

GOOCOOOOOOOOOOOOOOOOOOOGOOOOOOOOOOOOOOO

C

195

SUBRWTINE ISOTHERM
PARAMETER NK=5 ,NM=10,II=11.JJ=IO
CHARACTER‘ZO LABEL(NM)
COMMON IC21/XC(6.NK).YC(6.NK).ZC(6.6.NK)
IC22IDPHI.PHI(NM),ND(NK)
IC23IITSELF.IT,JT}KT.MT
IC24ITEM(II.JJ)
IC10IX(II).Y(JJ)
_ /C6/L.GTEM(8.250)
DATA DPHI.NDI1.E-10.2.2.4.2.6I
XC(1.1)=K(1)
XC(2.1)=X(2)
KC(1.2)=X(2)
XC(2.2)=X(3)
XC(1.4)=X(3)
XC(2.4)=X(5)
KC(1.5)=X(5)
KC(2.5)=X(6)
XC(3.5)=X(7)
KC(5.5)=K(10)
XC(6.5)=X(11)
DO 16 K=1.2
YC(1,K)=Y(1)
YC(2.K)=Y(2)
16 CONTINUE
YC(1.4)=Y(5)
YC(2.4)=Y(6)
YC(2.5)=Y(6)
YC(3.5)=Y(7)
YC(4.5)=Y(8)
YC(5.5)=Y(9)
DO 10 I=2.4
DO 10 J=1.4
XC(I.3)=X(I+3)
YC(J.3)=Y(J)
ZC(I.J.3)=TEM(I+3.J)
10 CONTINUE

XC(1,3)=X(3)
DO 20 1:114
ZC(1,J.3)=TEM(3.J)
20 CONTINUE
C
CW”INWT ZC DATA
C
no 7 I=1,2
DO 7 331,2

ZC(I.J.1)=TEM(I.J)

191+

ZC(I.J.2)=TEM(I+1.J)
ZC(I. J, 4)=T'EM(2‘I+1 .J+4)
7 CONTINUE
D0 8 I=4 . 6
D0 8 J=1 .6
ZC(I. J.5)*T'EM(I+5 .J+4)
8 CONTINUE
DO 9 I=1 .3
DG 9 J=1 .6
ZC(I.J.5)=TEM(I+4 .J'+4)
9 ENTINUE
C
Cm—DRAW ENTGUR OF COG. ING PLATE
C
CALL REGGVR
CALL G2DENT1( L)
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
C ICCCCCCCC IN ITIALIZATION BLOCK BLOCK 0100
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CPCCCCCCCC PROCESS BLOCK BLOCK 0200
C
C------DRAW IS G'IEERM 0F DI FFERET TEMPERATURE(K=1 . NK) AND DI FFERET
C----RE ION( M=1 . NM)
C
D0 1 K=1 , NK
DO 2 M=1 , NM
LAB=G
C
C----F IN D mE STARTING POINT OF IS OTHERM 0N UPPER BWNDARY
C
D0 4 I=1 .ND(K)-1
IF ((PHI(M)-ZC(I.ND(K).K))‘(PHI(M)-ZC(I+1.ND(K).K))
_ . LE.DPHI) TEEN
DX=XC (1+1 .K)-KC (I. K)
DZ=ZC(I+1 .ND(K) .K)-2N1. ND(K) .K)
XI=(PIII(M)-ZC(I. ND(K) .K) )‘DK/DZ+KC(I. K)
YI=YC(ND(K) .K)
CALL IDVEMXI . YI)
CALL SIDET( I . ND(K) . K. M)
12 IF (ITSEF.E.1) THE
CALL SIDEB( IT. 11'. KT. MI)
GO TO 12

ELSE IF (IT'SEF.E.2) TEE
CALL SIDET(IT. J’T.H.MT)
GO TO 12

195

ELSE IF (ITSEF.E.3) THEN
CALL SIDEL(IT.JT.KT.MT)
GO TO 12

ELSE IF (ITSEF.E.4) THE
CALL SIDER(IT.JT.KI'.MI)

 

GO TO 12
ED IF
LAB=1
GO TI) 17
END IF
4 CONTINUE
C
C—m-T-FIND THE STARTING MINT OF ISOTBERM ON THE RIGHT BWNDARY
C
17 DO 6 J=1.ND(K)-1
IF ((PRI(M)'ZC(ND(K).J.K))‘(PEI(M)“ZC(ND(K).J+1.K))
_ .LE.DPHI) THE
DY‘YC(J+19‘).YC(JOK)
YI=(PHI(M)-ZC(ND(K).J.K))‘DY/DZ+YC(J.K)
CALL IDVEA(XI.YI)
CALL SIDER(ND(K):J.K,M)
14 IF (ITSE.F.E.1) TEEN
CALL SIDEB(IT.JT.U.MT)
GO TO 14
ELSE IF (IT'SEF.BQ.2) THE
CALL SIDET(IT. TL“. '1')
GO TO 14
ELSE IF (ITSEF.m.3) THE
CALL SIDEL(IT.JT.KT.MI')
GO TO 14
ESE IF (ITSE.F.&.4) THE
CALL SIDER(IT.JT.KI.MT)
GO TO 14
ED IF
LAB=1
GO TO 18
ED IF
6 CONTINUE
C
Cum-FIND THE STARTING MINT OF ISOTHERM ON LWER BOUNDARY
C
18 DO 3 I=1.ND(K)"1
C
C-'-—-CHECK WHETHER TEE STARTmG HUNT OF ISOTHERM IS BETWEEN TWO
C- NODE
C

IF ((PEI(M)-ZC(I.1 .K) )‘(PHI(M)-ZC(I+1 .1 .K) ) .LE.DPHI)
THE

196

C

C-WCALCULATE COORDINATE OF STARTING NINT

C
DXBXC(I+1.K)’XC(I.K)
DZ=ZC(I+1.1.K)'ZC(I.1.‘)
XI‘(PHI(N)’ZC(I.1.K))‘DX/DZ‘FXCCLK)
CALL DVEA(XI.YI)
CALL SIDEB(I.1.K.N)

C

C----CHECK VALUE OF ITSEF AND DECIDE NEXT STEP
C
11 IF (ITSEF.m.1) mm
CALL SIDEB(IT.JT.U.IT)
GO TO 11
ESE IF unnama) mm
CALL SIDET(IT.JT.U.NI')
GO TO 11
ESE IF (ITSEF.m.3) THEN
CALL SIDE(IT.JT.H.NI')
GO '10 11
ESE IF (ITSEF.m.4) THE!
CALL SIDER(IT.JT.KT.UI')
GO TO 11
END IF
LAB=1
GO TO 19
MD IF
3 CONTINUE
C

C----FIND TEE STARTING POINT OF ISO'IBERM ON LEFI‘ BGJNDARY

C
19 DO 5 J=1.ND(K)-1

IF ( (PEI(N)‘ZC(1.J.K) )‘(PHI(M)-ZC(1 .J'+1 .K) ) .LE.DPHI)

HEN
XIBXC (1 .K)
DY=YC(J+1 .K)-YCLI. K)
DZBZC(1 91+]. .K)'ZC(1.J. K)

YI=(PHI(M)-ZC(1,J.K) )‘DY/DZ+YC(J, K)

CALL “NEARLY”
CALL SIDEL(1.J.K.N)
13 IF (ITSEF.m.1) THEN
CALL SIDEB(IT.J'T.KT. NT)
GO TO 13
ESE IF (ITSEF.m.2) TEE}
CALL SIDET(IT.JT.KT.UI')
GO TO 13
ESE IF (ITSEF.m.3) THEN
CALL SIDEL(IT.JT.H.IT)
GO TO 13

197

ESE IF (ITSEF.m.4) THE]
CALL SIDER(IT.JT.ET. '1‘)

GO TO 13
[N D IF
LAB=1
DID IF
5 CONTINUE
IF (LAB.m.1) GO TO 15
GO TO 2
C
C---LABE 'IBE TEMPERA'IURE OF IS O'IBERM
C
15 IF (X.LE.3) mm
CALL DVElH-O .1 .‘0 .7)
WRITE(LABE(N) .101) P3100
CALL CEAR'I'K(LABE(H) .0 .7)
101 FORMAT(F5 .0)
ESE
CALL HWER(-0 .1 .0 .25)
WRITE(LABE(H) .101) PENN)
CALL (BAR'IX(LABE(N) .0 .7 )
END IF
2 CONTINUE
1 CONTINUE
CALL AN IDDE
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CERCCCCCCC ERROR HES SAGES BLOCK 0900
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CF CCCCCCCC FORMAT STATEMENTS BLOCK 9000
C C
C C

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C

RE'IU RN

E“! D

198

CDCCCCC SUBRCIJTINE DESCRIPI' ION CCCCC
C C
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CV CCCCCCCC VARIABLE IDENTIFICATION CCCCCCCC
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CSCCCCCCC ENTRY _S'I’ORAGE BLOCK BLOCK 0000

C C

SUBROUTINE SIDET(I.J,K,N)
PARADETER NX=5.NM-10
COMLDN lC21/XC(6.NK).YC(6.NK).ZC(6.6.NK)

_ Iczz/npnl.m1(Nu).Nn(Nx)

_ Iczsnrsnm. IT, mm. II
c c
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c c
CICCCCCCCC INITIALIZATION BLOCK 13me 0100
c c
c c
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c c
cpcccccccc mocass BLOCK BLOCK 0200
c

nsmw=o

IF (1'1.LT.1) REJURN
IF ((PHI(H)“ZC(I+1.I-1.K))‘(PEI(N)’ZC(I+1.J.K)).LE.DPHI) TEEN
XI=XC(I+1.K)
DZ‘ZC(I+1,J.K)'ZC(I+1,I-1,K)
YI=(PEI(N)’ZC(I+1.I-1.K))‘DY/DZ+YC(J-1.K)
CALL RECOVR
CALL DRAWA(XI,YI)
CALL SIDEL(I+1,I-1,K.M)
ELSE IF ((PHI(N)-ZC(I.J-1.K))‘(PHI(N)’ZC(I+1.I’1.K)).LE.DPEI)
TEEN
DX=XC(I+1,K)-XC(I,X)
DZ‘ZC(I+1 o 1.1 ox)’ZC(I: 1.1 .K)
YI‘YC(I-1.K)
CALL RECOVR
CALL DRAWA(XI,YI)
I151
1151‘1
KT=K
"TB"

199

ITSELF=2
REIURN
ELSE IF ((PHI(N)-ZC(I,J-1.K))‘(PEI(N)-ZC(I,J.K)).LE.DPEI)
TEEN
XI=XC(I.K)
DY=YC(J.K)-YC(I-1.K)
DZ=ZC(I.J.K)-ZC(I.I-1.K)
YI=(PBI(N)-ZC(I.I-1.K))‘DY/DZ+YC(J-1,K)
CALL RECOVR
CALL DRAWA(XI.YI)
CALL SIDER(I.I-1,K,M)

END IF

c c
c c
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c c
canccccccc ERROR MESSAGES BLOCK 0900
c . c
c c
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c c
cpcccccccc FORMAT STATEMENTS BLOCK 9000
c c
c c

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C

REHURN

END

ZOO

CDCCCCC SUBRGJTINE DESCRIPTION CCCCC
C
C J(B : CONNECT ISO'IBERM FROM GIVE! POINT(CURSOR POSITION)
C TO THE POINT WHERE ISO'IEERN PASSING WE OF THREE O'IBER
C BCXJNDARY OF TH IS GRID .
C C
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CV CCCCCCCC VARIABLE IDENTIFICATION CCCCCCCC
C
C LJ :REFEREJCE NODE NUIBER OF THIS SUBRGJTINE
C
C K. M :REEION NUIBER AND ISO'lIIERM NUIBER WHICH 'IEIS SUBRGJTINE
C GOING '11) DRAW .
C
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CSCCCCCCC ENTRY _S'IORAGE BLOCK BLOCK 0000
C C

SUBROUTINE SIDEB( I . J. K. M)

PARAMETER NK=5 . NN=10

COMMN lC21/XC(6.NK) .YC(6.NK) .ZC(6.6.NK)

__ /CZZ/DPHI.PEI(NM).ND(NK)
_ IC23/ITSEF, IT, JT,KT,NI‘

C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CICCCCCCCC IN ITIALIZATION BLOCK BLOCK 0100
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CPCCCCCCCC PROCESS BLOCK BLOCK 0200
C

ITSEF=0
C

C---'-CHECK 'HE'IEER mE GRID IS WT OF RANGE OF IBIS REION
C
IF (J+1.GT.ND(K)) RE'IUW
C
C-WDRAW ISO'IEERM FROM GIVEN WINT TO THAT ON RIGHT BWNDARY
C
IF ((PHI(N)’ZC(I+1.J+1.K))‘(PHI(N)-ZC(I+1.J.K)).LE.DPHI) mm
XI=XC (1+1 . K)
DY=YC<J+1.K)*YC(J.K)
DZ=ZC(I+1,J+1.K)’ZC(I+1.J.K)

201

YI=(PBI(M)-ZC(I+1.J.K))‘DY/DZ+YC(J.K)
CALL RECOVR

CALL DRAWA(XI.YI)

CALL SIDEL(I+1.J.K.N)

C
C------ASSIGN A VALUE T0 ' ITSEF' IF ISO'IBERM PASSING UPPER
C----B(IJNDARY
C
ELSE IF ((PBI(l)-ZC( I. 1+1 .K) )‘ (PBI(N)-ZC(I+1 .J+1 .K) ) .LE.DPHI)
_ mm
DX=XC(I+1.K)-XC(I.K)
DZ‘ZC(I+1 0 1+1 3‘)-ZC(I: J+1 .K)
XII: (PHI(N)-'ZC(I. 1+1 .K) )‘DX/DZ-dC (I. K)
YIsYC(J+1.K)
CALL REGWR
CALL DRAWMXI . YI)
IT= I
JIbJ+1
KT=K
ITBM
ITSEF=1
mm
C
C----DRAW ISO'IBERM FROM GIVEN! POINT TO 'IBAT ON LEFT BWNDARY
C
ELSE IF ((PHI(M)-ZC(I, 1+1 .K) )‘ (PHI (M)-ZC(I, J. K) ) .LE. DPBI)
__ Tam
XI=XC (I . K)
DY=YC(J+1 .K)-YC (J. K)
DZ=ZC(I, 1+1 .K)-2N1. J'. K)
YI=(PHI(M)-ZC(I. J. K) )‘DY/DZ+YC(J. K)
CALL RECOVR
CALL DRAWMXI . YI)
CALL SIDER( I, J. K. M)
END IF
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CERCCCCCCC ERROR MESSAGES BLOCK 0900
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CFCCCCCCCC FORMAT STATEMENTS BLOCK 9000
C C
C C

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
RETURN

202

CDCCCCC SUBRWTINE DESCRIPI‘ ION CCCCC
C C
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CVCCCCCCCC VARIABLE IDENTIFICATION CCCCCCCC
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CSCCCCCCC ENTRY _STORAGE BLOCK BLOCK 0000

C C

SUBROUTINE SIDEL(I.J.K,N)
PARAMETER NK=5 . NIL-=10
COMMON IC21/XC(6.NK).YC(6.NK).ZC(6.6.NK)

IC22/DPEI,PBI(NN),ND(NK)
_ /C23/ITSELF.IT,JT,KT,NT
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CICCCCCCCC INITIALIZATION BLOCK BLOCK 0100
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CPCCCCCCCC PROCESS BLOCK BLOCK 0200
C
ITSELF=O

IF (I+1.GT.ND(K)) RE’IURN
IF ((PHI(M)-ZC(I,J.K))‘(PHI(M)-ZC(I+1,J.K)).LE.DPHI) THEN
DX=XC(I+1,K)-XC(I,K)
DZ=ZC(I+1.J,K)-ZC(I,J,X)
XI=(PHI(l)-ZC(I.J.K))‘DX/DZ+XC(I.K)
YI=YC(J.K)
CALL RECOVR
CALL DRAWMXLYI)
CALL SIDET(I.J.K.II)
ELSE IF ((PHI(M)-ZC(I+1.J.K))‘(PBI(M)-ZC(I+1.J+1.K)).LE.DPHI)
THIN
XI=XC (1+1 . K)
DY=YC(J+1,K)-YC(J,K)
=ZC(I+1.J+1,K)-ZC(I+1,J.K)
YI=(PHI(M)-ZC(I+1.J.K))‘DY/DZ+YC(J.K)
CALL RECOVR
CALL DRAWMXLYI)
ITbI+1
J’IbJ
KT=K
ur=u

203

ITSELF=3
RETURN
ELSE IF ((PHI(M)’ZC(I.J+1.K))‘(PRI(N)’ZC(I+1,J+1.K)).LE.DRHI)

THEN

DX=XC(I+1.K)’XC(I.K)
DZ=ZC(I+1pJ41.K)-ZC(I.J+1.K)
XI3(PHI(N)‘ZC(I,J+1,K))‘DX/DZ+XC(I.K)
YI=YC (1+1 9‘)

CALL RECOVR

CALL DRAIAIXI.YI)

CALL SIDEB(I.J+1.K.N)

END IF
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CERCCCCCCC ERROR MESSAGES . BLOCK 0900
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CFCCCCCCCC FORMAT STATEMENTS BLOCK 9000
C C
C C

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C

RESURN

END

204

CDCCCCC SUBRCIJTINE DESCRIPTION CCCCC
C C
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CVCCCCCCCC VARIABLE IDENTIFICATION CCCCCCCC
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CSCCCCCCC ENTRY _STORAGE BLOCK BLOCK 0000

C C

SUBRCIITINE SIDER( I , J. K. M)
PARADETER NK=5 , NMI=10
COMMON IC21/XC(6.NK).YC(6.NK).ZC(6.6.NK)

Iczzlnpn1.m1(Nu).Nn(Nx)
_ /c23/ ITSELF. IT. IT, IT, IT
c c
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c c
c1cccccccc INITIALIZATION BLocx BLOCK 0100
c c
c c
ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c c
cpcccccccc PROCESS BLOCK BLOCK 0200
c
ITSELF=0

IF (I-1.LT.1) RETURN
IF ((PHI(M)-ZC(I.J.K))‘(PEI(M)-ZC(I-1,J.K)).LE.DPHI) THEN
DX=XC(I.K)-XC(I-1.K)
DZ=ZC(I,J,K)-ZC(I-1,J,K)
XI=(PEI(M)-ZC(I-1,J.K))‘DX/DZ+XC(I-1,K)
YI=YC(J.K)
CALL RECOVR
CALL DRAMAIXI.YI)
CALL SIDET(I-1.J,K,M)
ELSE IF ((PEI(M)-ZC(I-1,J,K))‘(PEI(M)-ZC(I-1,J+1.K)).LE.DEHI)
TEEN
XIBXC(I-1.K)
DY=YC(J+1.K)’YC(J.K)
DZ=ZC(I-1,J+1,K)-ZC(I-1.J.K)
YI=(PHI(M)-ZC(I-1.J.K))‘DY/DZ+YC(J.K)
CALL RECOVR
CALL DRAIA(XI,YI)
IT%I-1
ITEJ
KT=K
MT=M

205

ITSELF=4
RETURN
ELSE IF ((PEI(M)-ZC(I.J+1.K))‘(PEI(M)-ZC(I-1.J+1.K)).LE.DPEI)

TEEN

DX=XC(I.K)-XC(I-1.K)
DZ=ZC(I.J+1.K)-ZC(I-1.J+1.K)
XI=(PEI(M)-ZC(I-1.J+1.K))‘DX/DZ+XC(I-1.K)
YI=YC(J+1.K)

CALL RECOVR

CALL DRAMAIXI.YI)

CALL SIDEB(I-1.J+1.K.M)

END IF
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CERCCCCCCC ERROR MESSAGES BLOCK 0900
C C
C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CFCCCCCCCC FORMAT STATEMENTS BLOCK 9000
C C
C C

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C

REHURN

END

206

CDCCCCC SUBRWTINE DESCRIPTION CCCCC

C

C 1(B :SUPPLY FORWARD DIFFEREICE mUATIONS 0F DIFFEREIT KINDS

C OF NODES SO THAT THEY CAN BE USE) TO CONSTRUCT THE

C WHCLE SET OF 'FDE' FOR A PARTICULAR SYSTEM.

C ('FDE': FORWARD DIFFERENCE NUATION)

C C

C C

C C

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC

C

CCCCCCCC VARIABLE IDBWTIFICATION CCCCCCCC

OOOOOOOOOOOOOOOOOOOOOOOGOOGOOOOOOGOGOOQO

RT:

DXI

DX2

DYl

DY2

II

12

11

12

TIME RATIO-TIN (Imm992/DIFFU SIVITY)

WHERE TDD IS TEE ARBITRARY TIME BASE ON WEIG TIME IS
NORMALIZED AND DIFFUSIVITY IS THAT OF THE MATERIAL OF
THIS NwE.

:NORMALIZED X-COORDINATE DIFFERBWE W.R.T. ITS RIGHT-'RAND
SIDE NODE.

:NORMALIZED X-COORDINAT'E DIFFERENCE. '.R.T. ITS LEFT-HAND
SIDE NODE.

:NORMALIZED Y-COORDINATE DIFFEREVCE M.R.T. ITS 'UPSTAIR'
NODE.

:NORMALIZED Y-COORDINATE DIFFEREQCE W.R.T. ITS
'DOWNSTAIR' NODE.

:X-DIRRECI'ION SUBSCRIPT OF NODE NUDBER FROM WHICH DO LOOP

IS STARTING .

:K-DIRRECI‘ION SUBSCRIPT OF NODE NUDBER AT WHICH DO LOOP

IS RIDING .

:Y-DIRRECTION SUBSCRIPT OF NODE NUDBER FROM WHIm DO LOOP

IS STARTING .

:Y-DIRRECI‘ION SUBSQIPT OF NODE NUIBH AT WHICH IX) LOOP

IS EQDING .

:CONDUCI‘ IV ITY RATIO=(COND. OF MATERIAL OF UPPER HALF”

(COND. OF MATERIAL OF LWER HALF) IN CASE OF COMNSITE
MATERIAL.

:DIFFUSIVITY RATIO=(DIFF. OF MATERIAL OF UPPER HALF”

(DIFF. OF MATERIAL OF LWER HALF) IN CASE OF COMPOSITE
MATERIAL.

207

C
C G :RATE OF HEAT GmERA'rION.
C
C UX :NORMALIZED AMBIENT TEMPERATURE ON LEFT on RIGHT SIDE.
C
C UY :NORMALIZED AIBIENT TEMPERATURE ON TOP OR BOTTOM SIDE.
C
C BIX :BIOT NULBER(=HEAT TRANSFER wEF.‘LmGIB/TBERMAL
C CONwCl‘IVITY) FOR LEFT OR RIGHT BCIJNDARY.
C
C BIY :BIOT NULBER<=HEAT TRANSFER COEF.‘LmG'm/THERMAL
C CONDUCTIVITY) FOR UP OR BOTTOM BwNDARI.
C
C C
c C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CSCCCCCCC ENTRY _STORAGE BLOCK BLOCK 0000
C C
C C
CcCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CICCCCCCCC INITIALIZATION BLOCK BLOCK 0100
C C
C c
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CPCCCCCCCC PROCESS BLOCK BLOCK 0200
C
C——---'FDE' FOR INNER NODES
C
SUBROUTINE mNER(RT,Dx1.Dx2,DY1,DY2.n.Iz,Jl,Jz,RK,RR,G)
COMIDN lCl/TRAX.TM.UN(11.10) .lCZ/DT.lu..UO(11.10)
Rx=Dx2/Dx1
RYsDn/Dn
DO 1 I=Il,12
DO 1 1811.12
OX=(UO(I-1.J)-(1.+1.IRX)’UO(I.J)+1.IRX‘UO(I+1,J))IDX1
IDx1
QY=(UO(I.J-1)-(1.+RK/RY)‘UO(I,J)+RK/RY‘UO(I.J+1))IDYI
IDYI
QG=G‘RL/ Dx1/ Dn
UN( 1.1) =2 .*DT*RT/ (1 .+RY¢RKIRR)‘( (1 .+RK‘RY) I (1 .+RX)‘QX+
_ RL‘RL‘QYH . «El (1 .+RK> )+UO(I. J)
1 CONTINUE
RETURN
END
C
C-----'FDE' FOR NODES 0N LEFT BOUNDARY

C

1

C

208

SUBROUTINE SLFI‘(RT. DX. DY1 ,DY2 . I, 11 .12 , RK, RR, G. UH,BIX)
COMDDN /C1/TMAX. TM. UN( 11 .10) ./C2/DT. IE. UO(11.10)
RY=DY2IDY1
D0 1 1811 .12
QX=(UO(I.1)-UO(I+1.1) )IDX/DX
QY=(U0(I. 1-1)’(1 .+RK/RY)‘U0(I. 1) +RK/RY‘UO(I. 1+1) ) IDY1/DY1
QG=G‘RL/DX/DY1
QM=(U0(I.1)-UH) IDX
UN( I . 1) =2 . ‘DT‘RT/ (1 .+RY‘RKI RR) ‘ (’( I .+RK‘RY) ‘QX+E'RL‘QY+
2 .‘QG—BIX'(1.+RY)‘QHX)+UO(I,1)

-CONTINUE

RETURN
ENID

C-----'FDE' FOR NODES (N RIGHT BGINDARY

C

1

C

SUBROUTINE SRGT(RT. DX, DYLDYZ . I, 11 .12 . RX, RR. G, UR.BIX)

COMM»! lC1/TMAX.TM.UN(11.10) .lC2/D'I‘.RL.UO(11.10)

RY=DY2IDY1

DO 1 J=Jl,12
QX=(UO(I-1.J)-UO(I.J))lDX/Dx
QY=(UO(I.J-1)-(1.+RK/RY)‘UO(I.J)+RK/RY‘UO(I.J+1))lDYI/DYl
QG=G‘RL/DX/DY1
QBX=(UO(I.J)-UR)/DX
UN( 1. .1) =2 .‘DT‘RT/ (1 .+RY¢Rx/RR)‘ ( (1 .+thRY)'0x+m.‘RL‘QY+

2.‘OG—BIX‘(1.+RY)‘OHX)+UO(I.J)

-CONTINUE

RETURN
BID

C------'FDE' FOR NODES m UPPER BWNDARY

C

1

C 7

SUBRWTINE SUP(RT,DX1 ,DX2 ,DY. I1 , I2 ,1, RR. RR, G. UH. BIY)
COMIDN lC1/TMAX.TM.UN(11.10) .lC:l/D'I‘.RL. UO(11.10)
RX=DX2IDXI
DO 1 I=I1,I2
QX=(UO(I‘1.J)‘(1.+RK/RX)‘UO(I.J) +RKIRX‘UO(I+1.J) )lDXl/DXI
QG=G‘RL/DX1/DY
QRY=(UO(I.J)’UH) IDY
UN( I . 1) =2 .‘DT‘RT/ (I J'n’RK/RR) ‘ (QX+(1 .+RK‘RX) ‘RL‘RL‘QY‘F
2 .QG-BIY.(I.+RX)‘RLQHY)+UO(I,J)

-CONTINUE

RETURN
END

C----'FDE' FOR NODES (N BOTTOM BCIINDARY

C

SUBRCIITINE SBTM(RT, DX1.DX2,DY, Il,IZ,1, RK,RR,G, UH, BIY)
COMDDN /C1/TMAX.TM.UN(11.10)./C2/DT.RL.UO(11.10)

209

RX=DX2IDX1

DO 1 1811.12
QY=(UO(I.J)'DO(I.J+1))IDY/DY
QX=(UO(I-1.J)'(1.+RK/RX)‘UO(I.1)+RK/RX’DO(I+1.J))IDn/Dn
QG=GPRLIDXIIDY
QHY:(UO(IOJ)“UB)IDY
UN( I . 1) =2 . .DT‘R'I/ (I .+n‘RK/RR) ‘ (OX-(1 .+RK‘RX)‘RL‘RL‘QY+

2.‘QG—BIY‘(1.+RX)‘RL‘QHY)+DO(I.1)

1 CONTINUE
RETURN
DID
C
C------ 'FDE' FOR NODE AT LEFT-DO'NIARD CORNER
C
SUBRCIITINE CLD(RT.DX,DY, I,J.G.UX.UY.BIX.BIY)
COMMON /C1/TMAI.TI.UN(11.10).IC2/DT.RL.UO(11.10)
QX=(UO(I.J)-U0(I+1.J))lDXIDx
QY=(UO(I.J)-UO(I.J+1))lDY/DY
QG=G‘RL/DX/DY
0HX=(UO(I.J)-UX)/Dx
QUE-(00(I.J)-UY)/DY
UN( I . 1) =2 . 'DT‘RT’ (~01-RL‘RL‘QY-BIX‘ORX-BIY'RL'QRY-t2 . ‘00)
1UO(I.J)
RETURN
END
C
C----- 'FDE' FOR NODE AT LEFT-WARD CORNER
C
SUBROUTINE CLU(RT. DX. DY. I. J. G. UX. UY. BIX. BIY)
COMLDN /C1/TMAX.TM, UN(11,10) .IC2/DT. RL. UO(11.10)
0x=(UO(I.J)-U0(I+1.J))IDXIDX
QY=(UO(I.J-1)-UO(I.1) )lDY/DY
QG=G‘RL/ DX/ DY
0HX=(UO(I.J)-UX)IDX
QHY=(U0(I.J)-UY)/DY
UN( I . 1) =2 . ‘DT‘RT‘ (-0.X+RL‘RL‘OY-BIX‘QRX-BIY‘RL‘QBY+2 . ‘00)
‘1U0(I.J)
RETURN
EJD
C
C-----'FDE' FOR NODE AT RIGHT-DOWNWARD CORNER
C

SUBRWTINE CRD(RT, DX. DY. I. 1. G. OX. DY. DIX. BIT)

COMN [CI/TEAL '1'“. UN( 11 .10) .lCZ/DT. RL. 00(11 .10)
QX‘(UO(I-1.J)-UO(I.J) )IDX/DX
QY'(UO(I.J)-UO(I.1+1))lDY/DY

QGBG‘RL/DXIDY

QHX=(UO(I.1)-UX)/DX

QHY=(UO(I.1)*UY)/DY

UN( I . J) ‘2 . .DT’RT‘ (QX’E‘RL‘QY‘BIX‘QHX'BIY‘RLQHY‘FE . .QG)

21C)

:U0(I,J)
RETURN
END
C
C-—--——'FDE' FOR NODE AT RIGHT-OPRARD CORNER
C
SUBROUTINE CRU(RT.DX. DY, 1.1. 0. UR. UY. BIX. BIY)
COMMON lC1/TMAX.TM.UN(11.10).lCZ/DT.RL.UO(11,10)
Qx=(UO(I-1.J)-UO(I.J))lDI/Dx
QY=(UO(I.J-1)-UO(I.J))lDY/DY
QG=G¢RLIDXIDY
QHX=(UO(I.J)-UX)/Dx
QEY==(UO(I.J)-UY)/DY
UN(I.J)=2.*DT¢RT*(0X+RLPRL‘QYEBIX‘QHx-BIYtnL‘QEY+2..QG)
1U0(I.J)
RETURN
END
C
C—-—-—---'FDE' FOR NODE AT LEFT-DORNWARD INNER CORNER
C
SUBRCIJTINE ICLD(RT. DX1 ,sz .DY1,DY2 , I, J, RR, RR, G, UX. UY. BIX. BIY)
COMMON /C1/TMAX.TM.UN(11.10).lC2/DT.RL.UO(11.10)
RX=DX2IDXI
RY=DY2IDY1
QI=(RK*RY*UO(I-1.J)-(RK¢RY+RK¢RY/Rx+1.lRX)‘UO(I.J)+(RK¢RY/Rx+
;./RX)‘UO(I+1.J))le1/Dx1
QY=(RX*UO(I.J-1)—(Rx+RK/RX+RK*RXIRX)’UO(I.J)+RK/R!*(1.+
RX) ‘UO(I . 1+1) ) IDY1/ DY1
QG=G‘RL/ DXl/ DYI
QRX=(UO(I.J)-UX)ID11
QHY=(UO(I.J)-UY)/DY1
UN( I . 1) =2 .‘DT’RT/ (RX+RY‘RK‘(1 .+RX) IRR) ‘ (OX+RL‘RL‘OY-BIX‘OEX-
gIYORL‘QHY+2.‘QG)+U0(I.J)
RETURN
END
C
C----—'FDE' FOR NODE AT LEFT-UNARD INNER CORNER
C
SUBRwTINE ICLU(RT,DX1.DX2.DY1,DY2,I.J,RK.RR,G,UX.UY,BIX,BIY)
COMMON lC1/TMAX.TM.UN(11.10).lC2/DT.RL.UO(11.10)
Rx=Dx2/Dx1
RY=DY2lDY1
Qx=(UO(I-1,J)-(1.+RK‘RY/Rx+1.IRX)‘UO(I,J)+(RKPRY/Rx+
;./RX)‘UO(I+1.J))IDx1/Dx1
OI=((1.+RX)‘UO(I.I—1)-(RX+1.+RK*RX/RI)‘UO(I.J)+RK/RI‘RX
gUO(I.J+1))/DY1/DY1
QG=GPRLI Dx1/ DY1
0nx=(UO(I.J)-UX)/Dx1
QEF=<UO(I.J)-UY)/DY1
UN(I.1)=2.‘DT'RT/(l.+RX+RI‘RK‘RX/RR)‘(0X+RL'RL‘QY‘BIX!OHX-

211

BIT-RL-Olmz . «10) +UO( I . J)
RETURN
mm
C
C-----—'FDE' FOR NODE AT RIGHT-UNABD INNER CORNER
C
SUBROUTINE ICRU(RT. Dn ,sz .011. DYZ . I. J, RK, RR, 0, U3. UY, BIx, BIY)
COMMON [Cl/TMAX.TM.UN(11.10) .IC2/DT.BL.UO(11.10)
Rx=Dx2/Dx1
RY=DY21DY1
(DH (1 .+RK‘RY)‘UO(I-1 .J)-(RKtRT+1 .+1 .ln)‘UO(I. J) +
1./ RX‘UO(I+1 .J) )/Dn/Dx1
QY=((1.+RX)’U0(I.1-1)-(RX+RKIRY+1.)*UO(I.J)+RKIRY‘
100(1 . J+1) ) Inn/011
00=0-RL/Dx1/Dn
ORI=(U0(I.J)-Ux)/Dx1
QBY=(UO(I.1)-UY)IDY1
um I . J) =2 . tDTtRT/ (1 .+n+RY‘RKl RR) . (mama-0141111111-
Bnan-mmz . ‘00) +U0( I. J)
RETURN
IND
C
C-----'FDE' FOR NODE AT RIGHT-DowNwARD INNER CORNER
c
SUBROUTINE ICRD(RT, 0x1 .sz .DY1 ,sz . I, J. RK. RR, 0. UL UY. BII, BIY)
comm lCl/TMAX.TM,UN(11,10) .IC2/DT.RL,U0(11.10)
RX=DX2IDX1
RY=DY2IDY1
OK=( (1 .+RK‘RY) ‘UO( I-1 , J) -(RKtRT+RKtRT/ Rx+1 . ) ‘UO(I , J) +RK¢RIIRK
gUO(I+1 .J) ) IDx1/Dx1
QY=(UO(I.J-1)-(1.+RK/RY+RK‘RXIRY)‘UO(I.J)+RKIRY‘(1.+
RX) ‘UO(I . J+1) ) /DY1/DY1
QG=G*RL/DXl/DY1
OBI=(UO(I.J)-Ux)/Dx1
QRY=(UO(I.1)-UY)IDY1
UN( I . J) =2 . ‘DT‘RTI (1 .+RK¢RK*(1 .+Rx) IRR) . (0X+RL'RL‘OY-BIX‘OEX-
BIY‘RL‘OHYH . .00) +U0( I. J)
RE'IURN
END
C c
C c
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C C
CERCCCCCCC ERROR MESSAGES BLOCK 0900
c C
C c
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C c
CFCCCCCCCC FORMAT STATEMENTS BLOCK 9000
C C

212

C C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C

RE’IU RN

END

10.

11.

REFERENCES

K. J. Henle and L. A. Dethlefsen. ”Time-Temperature Relationships
for Heat-Induced Killing of Mammalian Cells". in: ermal

Characteristics 2; Tumors: Anglications 1; Detection and
Treatments. ed. by R. K. Jain and P. I. Gullino. Ann. N.Y.
Acade 801.. ' v0103359 19809 ppe235-2350

J. J. MCGrath. The Dynamics of Freezing and Thawing Mammalian Cells:
the Hela Cell. I.S.I.E. Thesis. NIT. 1974

Vorapot Khompis. A.Numerical Heat Transfer Analysis of a CryomicroscOpe
Conduction Stage. H.S.H.E. Thesis. ISU. 1980

M. W. Scheiue and CH. Korber. Thermally Defined Cryomicroscopy and
Some Applications on Human Leucocytes. Journal of Microscopy.
Vol. 126. Pt 1. April 1982

A. Goldsmith and T. E. Waterman and H. J. Hirschhorn, Handbook of
ThermOphysical Properties of Solid Haterials. Pergmmon Press. 1962

J. P. Holman. Heat Transfer. Ichraw-Hill. Inc.. 1976

W. H. Rays and M. E. Crawford. Conyective Heat and Mass Transfer.
HCGraw-Hill. Inc.. 1966

Robert Jorgensen. Fan Engineering. Buffalo Forge Co.. 1961

G. J. Van Wylen and R. E. Sonntag. Fundamentals of Classical
Thermodynamics, John Wiley and Sons. Inc.. 1974

G. E. Myers. Analytical Methods in Conduction Heat Transfer,
McGraw-Hill. Inc., 1971

Omega Engineering Inc.. The Omega 1981 Tbmperature Measurement
Handbook. Omega Press. 1981