DINSRDC - 716/032 DAVID W. TAYLOR NAVAL SHIP RESEARCH AND DEVELOPMENT CENTER Bethesda, Md. 20084 DTNSRDC-78/032 A NONLINEAR MATHEMATICAL MODEL OF MOTIONS OF A PLANING BOAT IN REGULAR WAVES by Ernest E. Zarnick a e a WHOL ON DOCUMENT ) \ COLLECTION s, ° of APPROVED FOR PUBLIC RELEASE: DISTRIBUTION UNLIMITED SHIP PERFORMANCE DEPARTMENT RESEARCH AND DEVELOPMENT REPORT "MATHEMATICAL MODEL OF MOTIONS OF A PLANING BOAT IN REGULAR WAVES March 1978 DTNSRDC-78/032 MAJOR DTNSRDC ORGANIZATIONAL COMPONENTS DTNSRDC COMMANDER 0 TECHNICAL DIRECTOR 01 OFFICER-IN-CHARGE OFFICER-IN-CHARGE CARDEROCK ANNAPOLIS SYSTEMS DEVELOPMENT DEPARTMENT 1 SHIP PERFORMANCE te eee ee SURFACE EFFECTS a DEPARTMENT 46 COMPUTATION, STRUCTURES MATHEMATICS AND DEPARTMENT aa LOGISTICS DEPARTMENT 18 PROPULSION AND AUXILIARY SYSTEMS SHIP ACOUSTICS DEPARTMENT 47 DEPARTMENT 16 CENTRAL SHIP MATERIALS INSTRUMENTATION ENGINEERING DEPARTMENT 99 DEEAR MEN gs 0 CT 0 INTO IMMA | UNCLASSIFIED SECURITY CLASSIFICATION OF THIS PAGE (When Date Entered) a READ INSTRUCTIONS al | REPORT DOCUMENTATION PAGE | BEFORE COMPLETING FORM —_| . REPORT NUMBER 2. GOVT ACCESSION NO,| 3. RECIPIENT'S CATALOG NUMBER | 4. TITLE (and Subtitle) 5. TYPE OF REPORT & PERIOD COVERED A NONLINEAR MATHEMATICAL MODEL OF MOTIONS OF A PLANING BOAT 8. CONTRACT OR GRANT NUMBER(e) IN REGULAR WAVES » AUTHOR(a) Ernest E. Zarnick 10. PROGRAM ELEMENT, PROJECT, TASK AREA & WORK UNIT NUMBERS ZF 43 421001 Work Unit 1-1500-100 12, REPORT DATE March 1978 13. NUMBER OF PAGES 86 1S. SECURITY CLASS. (of thie report) UNCLASSIFIED 1Sa. DECL ASSIFICATION/ DOWNGRADING | SCHEDULE APPROVED FOR PUBLIC RELEASE: DISTRIBUTION UNLIMITED ” PERFORMING ORGANIZATION NAME AND ADDRESS David W. Taylor Naval Ship Research and Development Center Bethesda, Maryland 20084 - CONTROLLING OFFICE NAME AND ADDRESS Naval Sea Systems Command (SEA 035) Washington, D.C. 20362 - MONITORING AGENCY NAME & AODORESS(If different from Controlling Office) » DISTRIBUTION STATEMENT (of this Report) . DISTRIBUTION STATEMENT (of the abatract entered in Block 20, if different from Report) - SUPPLEMENTARY NOTES - KEY WORDS (Continue on reverse aide if neceseary and identify by block number) Planing Boat Motions Hydrodynamic Impact Small Boat Worthiness Nonlinear Ship Motions in Waves . ABSTRACT (Continue on reverse aide if neceseary and identify by block number) A nonlinear mathematical model has been formulated of a craft having a constant deadrise angle, planing in regular waves, using a modified low-aspect-ratio or strip theory. It was assumed that the wavelengths would be large in comparison to the craft length and that the wave slopes would be small. The coefficients in the equations of motion were determined by a combination of theoretical and empirical relationships. A simplified version for the case of a craft or model being towed at constant speed was programed for computations on a digital computer, and the results were compared with existing experimental data. (Continued on reverse side) DD en, 1473 EDITION OF 1 Nov 65 1s OBSOLETE 1 JAN 73 Pritciaate oilele sor UNCLASSIFIED SECURITY CLASSIFICATION OF THIS PAGE (When Data Entered) UNCLASSIFIED SECURITY CLASSIFICATION OF THIS PAGE (When Data Entered) (Block 20 continued) Comparison of computed pitch and heave motions and phase angles with corresponding experimental data was remarkably good. Comparison of bow and center of gravity vertical accelerations was fair to good. ; UNCLASSIFIED SECURITY CLASSIFICATION OF~THIS PAGE(When Data Entered) TABLE OF CONTENTS Page NIB SIRRIA Gli rasscts crs cionlatoar Ata S he a Sale oe Hes ELSE ac REE ae Ta a eS RAE tt ee eB MS hai ante tee ] ADMINISBRAGIVESINEORIMATIONP eae cecscsartese reer crecet eco eeacncccctaacstadnecbertasteeascae seater nes 1 INE RO DU GTO Neer cree ee asc hh esc sea oe ee ee a ea tae eet endataets Mebaadace UA taeen 2 SSAA ] MATHEMATICAT: FORMULATION 2.235 2R i Se emesenceietacaueendecctenetedseeteeuseteteuneds sbenken ct 2 GEINEIRIAT eels a Mic ome Miss kaka ry Riad eo eb da Mette Ay uattae dated ute aaa a lula cena u nae D TWO-DIMENSIONAL HYDRODYNAMIC FORCE ...........cceecccecceceeseetseseeseeeseeeeeeecees 3 TOTAL HYDRODYNAMI©G FORCE AND) MOMENT) «ore. ccc. ciccccccccescceescc-coccecececseeeee 6 EOQUATIONSIORSMOMON GENERA ee seere ere cee catea ee aeeraaes eonce eee eee een asaee eee 8 EQUATIONS OF MOTION, SIMPLIFIED FOR CONSTANT SREED iain Leen rt. 0 eae Peete sec nora bt vavalanana tneabisnaiestceenceveeteecenters 9 COMPARISON OF COMPUTED RESULTS WITH EXPERIMENTS ............ eens 10 CONCEUSIONSFAINDIRE COMMENDATION S eesccccceecrcsescoaseccesereceee cece cee ceeeeeccteeeceeeeeece 13 INCKNOWEEDG MENTS tocecio lta toy sesui Miaka nebete esta ccrurvumal Anse dbestanmccastn rat Sat hmalataernics Mlseaemvas 13 RO IEIREINGES riees ucuecasensl ati sausunasaauleiatsh Sas detest cdtewealel nasteleneatetaccuatien ar wrateoe dees meacetericmeatemstnc 333} APPENDIX A — EVALUATION OF HYDRODYNAMIC FORCE AND IMOMENTAINEGRAICS crete nce rerscecoacterauettcnedeadccussstmeceemeceu bean sachet 35 APPENDIX B — COMPUTER PROGRAM DESCRIPTIONS 1.0.0.0... ceeeeeeeceeeeeeeeeeeeeeeeeeeees Bo) LIST OF FIGURES a Coordinates SVstennescascccssrcesectsce desetescctesnesuemuascersboealst ences seelrele denamstetusweacecebauclsnesaieiomeistses 14 2, = TOES OI? TUNyrO> Diba aSOrmell JEON, bassssscoccqcccononssocosopsasanastonodbsadenbaocqosnaonansonsoadasedecoq0000 14 BE—WIbinesvOlmenSmMationMOdelSiawees.cayetsstantacauaaseccs usc Susees lawamcld somes sneauaeteswellesvesscamesuientsdess 15 4 — Sample Time Histories of Computed Pitch and Heave Motions...............cceeeeeeees 16 5 — Sample Time Histories of Computed Accelerations of Bow arrade Ge nite rao fit Gravity ieee ean ck creer ese ce eee a ratte Natrol aac DA Oe ee en eA 17 6 — Variation of Pitch and Heave with Wave Height 0.0... eeececccse esse eeeeeeeeeeneeeeee 18 7 — Variation of Acceleration of Bow and Center of Gravity AVIGHWAWIAV.CR El ele (ite ee a ouMRMe Noel NUE ee IN UES ae eal aati ary sav bualsusbinceece va usaubesesenuccildeceesser eas 19 ili 8'= Trajectory of Computer Model Relative to Wave ..............0.c.c.0..-ccessseserssooococernrocooens 20 9 — Heave Response for 10-Degree Deadrise Model at Wiki = Oe Ohne tetaanawcnica ete moet 21 10 — Pitch Response for 10-Degree Deadrise Model at V/./L= CLO eee: np 11 — Heave Response for 20-Degree Deadrise Model at Vik/ile ='6.0) coco ees 23 12 — Pitch Response for 20-Degree Deadrise Model at Wa = GiOb Net eet eee ee eee 24 13 — Heave Response for 30-Degree Deadrise Model at Wi ale = Oi in anes Ds 14 — Pitch Response for 30-Degree Deadrise Model at V/./L= O2Obscn cc aucamecsutace serene 26 15 — Heave Response for 20-Degree Deadrise Model at Min) = 4 Ole tea ee a a OF) 16 — Pitch Response for 20-Degree Deadrise Model at VIS L S40 FNL OME oe eee 28 17 — Bow Acceleration for 10-Degree Deadrise Model at Wy b= 610 eee ee 29 18 — Center of Gravity Acceleration for 10-Degree Deadrise Model FIL 0) 0 Ow tO ninarereener a Ra Nerealcate eM edt dE ae adc ae etch i sxotcscqoo0 30 19 — Bow and Center of Gravity Accelerations for 20-Degree Deadrise) Modelvatwvi/y/ale— 4+0) amdVi//x/ ole. — (61 Ole ercetes enone eee eee 31 20 — Bow and Center of Gravity Accelerations for 30-Degree DeadriseModeliatyVin/Ml OO ce atu eee cle ent 1 ee ke 32 Table 1 — Model Characteristics and Wave Conditions for Computations ......................55 11 NOTATION Mass matrix Section area Correction factor for buoyancy force Half-beam of craft Crossflow drag coefficient Load coefficient A/pg(2b)3 Wavelength coefficient L/A [Ca /( 2b)2 v2 Friction drag force Total hydrodynamic force in x direction Total hydrodynamic force in z direction Total hydrodynamic moment about pitch axis Two-dimensional hydrodynamic force Acceleration of gravity Wave height, crest to trough Vertical submergence of point below free surface Double amplitude of heave Pitch moment of inertia Added pitch, moment of inertia Wave number Two-dimensional added-mass coefficient Hull length Longitudinal center of gravity, percent of L Mass of craft Added mass of craft ~~ a Sectional (two-dimensional) added mass Hydrodynamic force normal to baseline Wave elevation r= Ty COS (kx+wt) Wave amplitude Relative fluid velocity parallel to baseline Relative fluid velocity normal to baseline Speed-to-length ratio in knots/ft!/? Weight of craft Vertical component of wave orbital velocity Vertical component of wave orbital acceleration Fixed horizontal coordinate Vector of state variables Surge velocity Surge acceleration Surge displacement Fixed vertical coordinate Heave velocity Heave acceleration Heave displacement Deadrise angle Hull displacement W Body coordinate normal to baseline Wavelength Pitch angle Pitch angular velocity vi Pitch angular acceleration Double amplitude of pitch Body coordinate parallel to baseline Density of water Wave frequency Wetted length Vii 4 ay at 4 Des bade , hi, i } Hy) os ’ , a i 1 Na ‘ay t indy iit, Rapes a ‘ ey. P > Wises, ahavade we Wad Bice Baa t F ; Wiha sei by oe Ki Fs < = a Like Rave "CG et) i & ra * ABSTRACT A nonlinear mathematical model has been formulated of a craft having a constant deadrise angle, planing in regular waves, using a modified low-aspect-ratio or strip theory. It was assumed that the wavelengths would be large in comparison to the craft length and that the wave slopes would be small. The coefficients in the equations of motion were determined by a combination of theoretical and empirical relationships. A simplified version for the case of a craft or model being towed at constant speed was programed for computations on a digital computer, and the results were compared with existing experi- mental data. Comparison of computed pitch and heave motions and phase angles with corresponding experimental data was remarkably good. Comparison of bow and center of gravity vertical accelerations was fair to good. ADMINISTRATIVE INFORMATION This investigation was authorized by the Naval Sea Systems Command with initial funding under Task Area SR-023-0101 and completion under Task Area ZF-43-421001. INTRODUCTION Computer programs for estimating the motions of displacement ships in waves for all headings and speeds have been in existence for some time. Comparable computational schemes for planing craft do not exist except in limited and restricted cases. A program for planing craft would be quite useful to the small craft designer, providing a means for systematically exploring the effects of numerous design variations on performance of the craft in waves. With minor modification, the program could also be used to examine the merits of a hybrid craft design, e.g., a combination of planing craft and hydrofoil. Predicting the motions of a planing craft in wave’s is by no means a simple problem. The analytical description of a high-speed craft, planing in waves, involves several different types of flow phenomena, including planing; hydrodynamic impact, and, to a lesser extent, surface wave generation and hydrostatics. Also, the mathematics tend to become nonlinear rapidly as the motion increases or, like the real craft, can in some instances exhibit large instabilities such as porpoising. Development of a computer program that would take into account all of the previously described factors and would be applicable for a wide range of speed and wave conditions requires a careful and systematic study in several stages with appropriate verification at each stage. To lay the foundation for such a general program, a simpler problem has been formulated in this report with potential for expansion and generalization to the more complicated case. The simpler problem is that of a V-shaped prismatic body with hard chines and constant deadrise planing at high speed in regular head waves. The mathematical formulation is analogous to low-aspect-ratio wing theory with provisions for including hydrodynamic impact loads, essentially a strip theory. Surface wave generation and forces associated with unsteady circulatory flow are neglected, and the flow is treated as quasi-steady. The mathematical formulation is an empirical synthesis of several theoretically derived flows describing the overall craft hydrodynamics. Wave input is restricted to monochromatic linear deepwater waves with moderate wavelengths and low wave slopes. MATHEMATICAL FORMULATION GENERAL Consider a fixed coordinate system (x,z) (Figure 1) with x axis in the undisturbed free surface, pointing in the direction of craft travel, and the z axis, pointing downward. If the motions of the craft are restricted to pitch 0, heave Zcq> and surge Xcq, the equation of motions can be written as MXcg = T,,-N sin @.—Djcos 0 MZcc = T,-Ncos#+Dsin@ +W 10 = Nx, + Dx Tx, (1) where M is mass of craft I is pitch moment of inertia of craft N_ is hydrodynamic normal force D_ is friction drag W is weight of craft is thrust component in x direction is thrust component in z direction x, is distance from center of gravity (CG) to center of pressure for normal force Xq is distance from CG to center of action for friction drag force x. is moment arm of thrust about CG. Equation (1) is exact; however, defining the hydrodynamic forces and moments in waves can be extremely difficult. A high-speed craft moving in waves may transit through several regimes that have different hydrodynamic flow characteristics. For example, as the craft moves away from the crest of wave, the flow may be characterized by unsteady-state planing until the craft collides with the oncoming wave crest and enters another regime in which impact forces are important. After the impact, the craft may enter still another regime in which it is planing but in which buoyancy forces are rather significant. The most promising approach to a method that would incorporate all three types of flow conditions into a general formulation would seem to be a modified strip theory. The mathematical justification for this approach is not rigorous; however, there is sufficient precedent to expect promising results. For example, impact loads on landing seaplanes can be estimated reasonably well using a strip theory incorporating the Wagner two-dimensional (2-D), expanding-wedge theory,! and Chuang? has provided a strip method for determining loads on an impacting prismatic form that agrees extremely well with experimental results. More recently, Martin? has developed a linear strip theory for estimating motions of a planing craft at high speed, which shows good agreement with experimental results. A nonlinear model of the equations of motion would be expected to provide, in addition to the motions, reasonable estimates of the vertical accelerations which are an important consideration in designing a planing craft. TWO-DIMENSIONAL HYDRODYNAMIC FORCE Implicit with any strip method is the need to define the 2-D hydrodynamic force acting on an arbitrary cross section of the body. The 2-D flow problem is not simple; however, it lends itself to an empirical approach, using a combination of techniques used in hydrodynamic impact and low-aspect-ratio theories. The typical cross section of a hard-chine, V-shaped prismatic body such as that being considered here is shown in Figure 2. Figure 2 actually illustrates two different idealized- flow conditions, assumed to represent the crossflow during unsteady planing, depending upon whether the flow separates from the chine (Figure 2a) or not (Figure 2b). Nonwetted-chine flow conditions are typical of the sections near the leading edge of the wetted length of the craft. Wetted-chine flow conditions are more typical of sections near the stern, except possibly in the most extreme motion and wave conditions. Some sections between leading edge and stern may alternate between flow conditions as the wetted length changes with the motions. *A complete listing of references is given on page 33. 3 The normal hydrodynamic force per unit length f, acting at a section, is treated as quasi-steady and is assumed to contain components proportional to the rate of change of momentum and the velocity squared (drag term), i.e. Pes all Gate: 2 =~ |p (ma) + Ope bVv (2) where V is the velocity in plane of the cross section normal to the baseline m, is the added mass associated with the section form C,, . is the crossflow drag coefficient D,c p is the density of the fluid b is the half beam. For sections near the leading edge of the wetted length with nonwetted chine, the added mass is assumed to be defined in the same manner as during an impact which for a V-shaped wedge is given by m, =k, 7/2 pb? (3) where k, is an added-mass coefficient that may also include a correction for water pileup- k, is assumed to be 1.0 without pileup correction. The rate of change of momentum of the fluid at a section is given by io) dg dé D alee ay DS Pid eae (Ca De (4) where & is the body coordinate parallel to the baseline; see Figure 1. The last term on the right-hand side of Equation (4) takes into account the variation of the section added mass along the hull. This contribution can be visualized by considering the 2-D flow plane as a substantive surface moving past the body with velocity U = -dé/dt tangent to the baseline. As the surface moves past the body, the section geometry in the moving surface may change with a resultant change in added mass. This term exists even in steady-state conditions and is the lift-producing factor in low-aspect-ratio theory. The added mass of a section with fully wetted chines has not been developed to the same extent as the V wedge. In steady-state planing problems such as those of Shuford,* the crossflow is treated as a Helmholtz-type flow in which the Bobyleff results are used for estimating drag coefficients. Helmholtz flows are applicable only to steady-state conditions; so, it is assumed that the added mass for the fully wetted chine flow can be determined from Equation (3) using the value of the half-beam at the chine. In using the Shuford approach, it is assumed that the crossflow drag coefficient for a V-section is equal to the drag of a flat plate (C, , = 1.0) corrected by the Bobyleff flow coefficient approximated by cos 8, i.e. Coc = 1.0 cos 8 (5) The Bobyleff flow coefficient is the theoretical ratio of the pressure on a V-section to that experienced by a flat plate for a Helmholtz-type flow. The same approximation is used for estimating the drag coefficient for nonwetted chine sections, using the instantaneous value of the half-beam at the free surface. An additional force acting on the body is the buoyancy force f,,- This force is assumed herein to act in the vertical direction and to be equal to the equivalent static buoyancy force multiplied by a correction factor, i.e. fi =tee(A) (6) where A is the cross-sectional area of the section, and a is a correction factor. The full amount of the static buoyancy is not realized because at planing speeds the water separates from the transom and chines, reducing the pressure at these locations to atmospheric or less than the equivalent hydrostatic pressure. A greater reduction is realized in the buoyancy moment because of the corresponding shift in the center of pressure. Shuford4 in his work on steady-state planing recommended a factor of one-half to obtain the correct buoyancy force. In the following computations, the buoyancy force was corrected by a factor of one-half, ie., a= 1/2. The buoyancy moment, computed as the static buoyancy force multiplied by its corresponding moment arm, was corrected by an additional factor of one-half to obtain the proper mean-trim angles. Equation (2) is a synthesis of several idealized flow conditions combined in an empirical manner. In all of these flows, it is assumed that the net relative movement of the fluid past the body is in an upward direction. This condition may not always be met in the case of unsteady planing in waves. Closer scrutiny will be required to determine what limitations will be imposed upon the problem as formulated and/or what modifications will be required to improve the formulation. TOTAL HYDRODYNAMIC FORCE AND MOMENT The total normal hydrodynamic force acting on the body is obtained by integrating the stripwise, 2-D, hydrodynamic force given by Equations (2) and (6) over the wetted length 2 of the body. A body coordinate system (£,¢) with its origin at CG and the & axis pointing forward parallel to the baseline of the body is defined in Figure | to facilitate this integration. The hydrodynamic force acting in the vertical or z direction of the fixed integral coordinate system is given by -Ncos@) = FU(t) = / f cos @dé +f igat Q Q = | ftraeoveo + m, (,t) V(E, t) Q a0} dg +Cy Deb EN VE,D} cos dé +apsade] (7) where the integration is taken over the instantaneous wetted length. Similarly the force 1 acting in the horizontal or x direction is given by F, = [tin oat = — [ {mge.0 ven + m,(é,t) V(t) - UE) = [m,(E.0 VED) + Cp .é.teb,0V7(E,b} sin 0 dé (8) Wave forces are obtained by neglecting diffraction and assuming that the wave excitation is caused both by the geometrical properties of the wave, altering the wetted length and draft of the craft, and by the vertical component of the wave orbital velocity at the surface w,, altering the normal velocity V. The horizontal component of orbital velocity is neglected, since it is assumed small in comparison with the forward speed arae The velocities U and V may then be written as U= Xog Cos 6 - (Z¢g-Wz) sin 0 V = Xcq sin @ - bg + (Z6g-Wz) cos 9 The depth of submergence h of the body at any point P(£,¢) may be determined by h = Zceg - & sin @ + § cos @ -r where r is the instantaneous value of the wave elevation directly above the point. For regular head waves the wave elevation for a linear deepwater wave is TS atcOs k(x+ct) where r, is the wave amplitude k is the wave number c is the wave celerity. At point P(é,¢) X = Xo, t Ecos 8 +§ sin 8 (9) (10) (11) (12) The hydrodynamic moment Fg about CG is obtained in a similar manner by integrating over the wetted length the product of the normal force per unit length and the corresponding moment arm. Fees - [re nece =i cos 6 &dé Q Q =f {me001E0 + ieOVED — ULE 3 (mE OVE) + Gy EDP bE DVED +apgA cos 6} Edé (13) EQUATIONS OF MOTION, GENERAL Integrating the first term in Equations (7), (8), and (13) provides hydrodynamic forces and moments proportional to acceleration of the motion. These can be combined with the inertial terms of the rigid body to give the following equation of motion (M+M, sin? 9) Xa +(M, sin @ cos Nae =(QFysin 0)0 = 1, +F, -Dcosd (14) (M, sin 6 cos 9) Xa +(M+M, cos” Ware -(Q, cos 0)0 =17,+F,+Dsin@+W -(Q, sin 0)X og -(Q, C08 M%og + (1+1,)4 = Fg -Dxg + Tx, where M,(t) = [ me.nay Q Q,(t) = of m,(£,t)Edé Q L(t) = We m, (£,t) £7 dé Q F, = F, -{-(M, sin? )%¢q -(M, sin @ cos 0)%og + (Q, sin 0)6} Ee alse {appropriate acceleration terms} Fg = Fo - {appropriate acceleration terms}. A detailed evaluation of the integral expressions for the hydrodynamic forces and moments is provided in Appendix A. The solution to Equation (14) is cumbersome; however, it can be accomplished using standard numerical techniques. Introducing the state vector [x;, X4,X3,X4. Xe, X¢] where X) = Voc Xy = 20g or ACG CELA X¢ = Equation (14) can be rewritten, using matrix algebra, as > > Ax = g (15) so that Te eA eal ee x=A’g (16) where A“! is inverse of the inertial matrix A. Equation (16) is now in a form that lends itself to integration by using a numerical method such as the Runge-Kutta-Merson integration routine. EQUATIONS OF MOTION, SIMPLIFIED FOR CONSTANT SPEED Assuming that the perturbation velocities in the forward direction are small in comparison to the speed of the craft, the equations of motion may be further simplified by neglecting the perturbations and setting the forward velocity equal to a constant, i.e. og = CONSTANT If it is also assumed that the thrust and drag forces are small in comparison to the hydrody- namic forces and that they are acting through the center of gravity, the equations of motion may be written as Kats iO (M+M, cos? 0)%., -(Q, cos0)0 = F, +W -(Q, C08 O)Zcg + (I+1,)6 = Fo These equations also represent the case of the craft (model) being towed through CG at CONSTANT speed. Based upon the previously described equations of motion, a computer program has been written in FORTRAN language to compute the motions of a prismatic body, planing in regular head waves at high speed. A listing of the program along with the appropriate flow chart is presented in Appendix B. The listing contains reference to thrust and drag terms; however, they have no significance, except to provide a starting point for possible updating of the program to include these terms in the future. COMPARISON OF COMPUTED RESULTS WITH EXPERIMENTS Computations of pitch and heave motions and heave and bow accelerations were made, using the computer program for comparison with the experimental results of Fridsma.> Fridsma tested a series of constant-deadrise models of various lengths in regular waves to define the effects of deadrise, trim, loading, speed, length-to-beam ratio and wave proportions on the added resistance, heave and pitch motions, and impact accelerations at the bow and center of gravity. Figure 3 shows the lines of the prismatic models. The models were towed at CG with a system that permitted freedom in surge. The computer program simulates the model being towed at constant speed with CG at the baseline. Table 1 presents some characteristics of the model and experimental conditions for which comparisons were made. Most of the comparisons have been made at a speed-to-length ratio Wile of 6.0 where the mathematical model is expected to be most representative. A limited comparison has also been made at Wifes = 4.0; however, no comparison has been made at Wilx/ ale = 2.0. At this speed, the model (or craft) operates in the displacement mode for which the mathematical formulation is not valid. The average computer run corresponded to 10-second, real-time, model scale; however, only the last 2 seconds were considered free of transient effects. An example of the compu- ter time histories of pitch and heave motions is shown in Figure 4. Although the motions are periodic, they are not perfectly sinusoidal; consequently, in determining phase relationship, the peak, positive-pitch value (bow up) and the peak, negative-heave value (maximum upward position of CG) were used as reference points. There was a difference when the opposite peaks were used. TABLE 1 — MODEL CHARACTERISTICS AND WAVE CONDITIONS FOR COMPUTATIONS (Model Length = 114.3 cm (3.75 ft); L/b = 5; Ca = 0.608) CONFIGURATIONS B LCG Radius of SYMBOL Gyration V/V L 40 20 20 10 30 WAVE CONDITIONS FOR CONFIGURATION —-— ee | eee Corresponding time histories of bow and CG accelerations are shown in Figure 5. The bow acceleration was computed at Station 0. As can be seen in these plots, the impact accelerations ranged in magnitude from cycle to cycle. The maximum impact (or negative value) acceleration computed during the final 2 seconds of run was used in the comparisons with experimental values. In some instances, particularly near resonance, the maximum impact acceleration was more than twice the average impact value. Figure 6 shows a comparison of variation of computed and experimental pitch and heave motion with wave height for the 20-degree deadrise model in a 15-foot wavelength and for a speed-to-length ratio of 6.0. Figure 7 shows the corresponding impact acceleration at the bow and CG. The computed results closely follow the experimental data, except for CG acceleration at the extreme wave height condition, where the computed value is apparently much lower. Experimental data show that the model was leaving the water at this wave- height condition. The computer model did not leave the water but came very close; 11 see Figure 8. Figure 8 is a trajectory of the computer model relative to the wave for a selected cycle of motion. The computer model behaves very much as expected. On the left- hand side of the figure, the craft is planing down the crest of the wave and, as it approaches the wave trough, comes very close to leaving the water before slamming and submerging itself deeply into the front of the oncoming wave crest. Figures 9 through 14 show comparisons of the computed and experimental pitch and heave motions at VIS L = 6.0 through a range of wavelengths and at a constant wave height of 2.54 centimeters (1 inch) for deadrise models with 10, 20, and 30 degrees. The data have been plotted with respect to the coefficient C), defined by Fridsma as L/A [Ca/(L/2b)7] M3) Note that in our notation, b is the half-beam. Comparisons of heave and pitch for the 10-degree deadrise model shown in Figures 9 and 10, respectively, show excellent results. The computer model accurately predicts the secondary peaks in the pitch and heave responses at C, = 0.19. At this condition, the physical experimental model rebounds so as to fly over alternate waves. The computer model oscillates at half the wave-encounter frequency and comes close to leaving the water at alternate encounters with the wave. It does not quite leave the water to fly over alternate wave crests; nonetheless, it is a good representation of the actual motion. The heave and pitch comparison for the 20-degree deadrise model at WH TE = 6.0 is also excellent as can be seen in Figures 11 and 12, respectively. No experimental phase data for the condition were reported for C) greater than 0.072; however, extrapolated results (not shown) are in line with the computed results. The pitch and heave results shown in Figures 13 and 14 for the 30-degree deadrise model are good; however, responses at C) = 0.048 and C) = 0.072 are higher than the experimental results. For practical considerations a computational scheme for planing boat motions should be valid for a range from approximately Wale = 4.0 to Vis = 6.0. Computations of the motions were made for W/E = 4.0 for the 20-degree deadrise model; see Figures 15 and 16. Again the comparison of the computed heave and pitch response with experimental results is excellent. Comparisons of the computed and experimental impact accelerations (or largest negative values) are presented in Figures 17 through 20. Figures 17 and 18 show bow and CG accelerations for the 10-degree deadrise model; Figure 19 shows similar results for the 20- degree deadrise model; Figure 20 shows the results for the 30-degree deadrise model. In all cases, the comparison appears to be fair to good. In the shorter wavelengths, A/L = 1.0 and A/L = 1.5, the computed accelerations are higher than the corresponding experimental values. This is most pronounced for the 10-degree deadrise angle model. 12 a CONCLUSIONS AND RECOMMENDATIONS A mathematical model of a craft having a constant deadrise angle, planing in regular waves, has been formulated using a modified low-aspect-ratio or strip theory. It was assumed that the wavelengths were long in comparison to the craft length and that the wave slopes were small. The coefficients in the equations of motion were determined by a combination of theoretical and empirical relationships. A simplified version for the case of a craft or model being towed at constant speed was programed for computations on a digital computer, and the results were compared with existing experimental data. The comparison of the computed pitch and heave motions and phase angles with the corresponding experimental data gave remarkably satisfying results. Comparison of the bow and CG accelerations was fair to good. In summary, the previously described mathematical model appears to be a valid represen- tation of a planing craft in waves for the specific craft geometry and wave conditions considered. To make the computer program more valuable to the designer the following additional work is recommended: 1. Improve estimates of hydrodynamic coefficients to obtain better acceleration data and to include more complicated ship geometry. 2. Determine added resistance in waves. 3. Include freedom to surge and to add components of propulsion. 4. Extend to the case of irregular waves. ACKNOWLEDGMENTS Acknowledgment is given to Dr. Joseph Whalen and Ms. Sue Fowler of Operations Research, Inc., who translated the equations of motion into an operational computer program. Figure 1 — Coordinate System Figure 2a — Flow Separation from Chine Figure 2b — Nonwetted Chine Figure 2 — Types of Two-Dimensional Flow 14 SNOILVLS L “~L sp A x WYOANV Id (¢ auaJaJoy WO) S]OPO| OCUISLIg JO SOUT 0G "G0 ',,9€ — € ons Z‘AL‘L ‘ze, 'Y, o0E 002 oOL 1S PITCH (deg) HEAVE (cm) PITCH (rad) 0.20 0.15 0.10 0.05 0 -0.05 B = 20° V/\/L = 6.0 0.30 A/L = 4.0 HEAVE (ft) r, = 5.08 cm (2 in.) 0.10 TIME Figure 4 — Sample Time Histories of Computed Pitch and Heave Motions 16 CG ACCELERATION (G’‘s) BOW ACCELERATION (G's) B = 20° V/,/L = 6.0 A/L = 4.0 H = 5.08 cm (2 in.) TIME Figure 5 — Sample Time Histories of Computed Accelerations of Bow and Center of Gravity 17 DOUBLE AMPLITUDE OF HEAVE (cm) DOUBLE AMPLITUDE PITCH MOTION (deg) (inches) LEAVING WATER A EXPERIMENTAL (REFERENCE 5) O COMPUTED 0.1 0.2 0.3 WAVE HEIGHT/BEAM Figure 6 — Variation of Pitch and Heave with Wave Height CG ACCELERATIONS (G’s) BOW ACCELERATIONS (G's) LEAVING WATER A EXPERIMENTAL (REFERENCE 5) O COMPUTED Vi/L=6 A/L=4 B = 20° 0.1 0.2 0.3 WAVE HEIGHT/BEAM Figure 7 — Variation of Acceleration of Bow and Center of Gravity with Wave Height SABM O} DANLTIY [PPO JojndwoD jo A1ojsalery — g andy HONOYL JAVM LS3YD JAVM ee v & PHASE LAG HEAVE (h,/H) 100 2.0 = oi oy So 0.5 0 0.05 0.10 A EXPERIMENTAL (REFERENCE 5) O COMPUTED 0.15 Cy 0.20 EXPERIMENTAL MODEL REBOUNDED SO AS TO FLY OVER ALTERNATE WAVES 0.25 Figure 9 — Heave Response for 10-Degree Deadrise Model at V/,/ L= 6.0 All 0.30 PHASE LEAD LAG A EXPERIMENTAL (REFERENCE 5) © COMPUTED PITCH (6 /2mH/N) 0 0.05 0.10 0.15 0.20 0.25 0.30 Cy Figure 10 — Pitch Response for 10-Degree Deadrise Model at V/,/ L = 6.0 PHASE LAG HEAVE (h,/H) A EXPERIMENTAL (REFERENCE 5) O COMPUTED = (=) 0 0.05 0.10 0.15 0.20 0.25 Cy Figure 11 — Heave Response for 20-Degree Deadrise Model at V//L= 6.0 23 PHASE PITCH (0, /2mr/A) 100 A EXPERIMENTAL (REFERENCE 5) O COMPUTED 2.0 1.5 1.0 0.5 0 0.05 0.10 0.15 0.20 0.25 Cy Figure 12 — Pitch Response for 20-Degree Deadrise Model at V/,/ L = 6.0 24 PHASE LAG A EXPERIMENTAL (REFERENCE 5) © COMPUTED z HEAVE (h_/H) ron) 0 0.05 0.10 0.15 0.20 0.25 Cy Figure 13 — Heave Response for 30-Degree Deadrise Model at V/,/ L = 6.0 25 A EXPERIMENTAL (REFERENCE 5) O COMPUTED PITCH (0, /2mH/A) 0 0.05 0.10 0.15 0.20 0.25 Cy Figure 14 — Pitch Response for 30-Degree Deadrise Model at V/,/ L= 6.0 26 PHASE LAG z HEAVE (h_/H) A EXPERIMENTAL (REFERENCE 5) © COMPUTED = (=) 0 0.05 0.10 0.15 0.20 0.25 Cy Figure 15 — Heave Response for 20-Degree Deadrise Model at V/,/ L = 4.0 PITCH (0_/27H/)) p A EXPERIMENTAL (REFERENCE 5) O COMPUTED 0 0.05 0.10 0.15 0.20 0.25 0.30 Cy Figure 16 — Pitch Response for 20-Degree Deadrise Model at V/,/ L= 4.0 28 BOW ACCELERATIONS (G’s) A EXPERIMENTAL (REFERENCE 5) O COMPUTED 0 0.05 0.10 0.15 0.20 0.25 0.30 Cy Figure 17 — Bow Acceleration for 10-Degree Deadrise Model at V/,/ L = 6.0 29 CG ACCELERATIONS (G’s) A EXPERIMENTAL (REFERENCE 5) O COMPUTED 0.05 0.10 0.15 0.20 0.25 Cy Figure 18 — Center of Gravity Acceleration for 10-Degree Deadrise Model at V/,/ L = 6.0 30 0.30 CG ACCELERATIONS (G’s) BOW ACCELERATIONS (G’s) A EXPERIMENTAL (REFERENCE 5) O COMPUTED 0.05 0.10 0.15 0.20 0.25 0.30 Cy Figure 19 — Bow and Center of Gravity Accelerations for 20-Degree Deadrise Model at V/,/ L = 4.0 and V/,/ L= 6.0 Sil CG ACCELERATIONS (G's) BOW ACCELERATIONS (G’‘s) A EXPERIMENTAL (REFERENCE 5) O COMPUTED 0.05 0.10 0.15 0.20 0.25 Cy Figure 20 — Bow and Center of Gravity Accelerations for 30-Degree Deadrise Model at V/,/ L= 6.0 32 REFERENCES 1. Wagner, H., “Landing of Seaplanes,”’ leitschrift fur Flegtechnik und Motorluflsshiff- fahrt, (14 Jan 1931); National Advisory Committee for Aeronautics TM 672 (May 1931). 2. Chuang, S.L., “Slamming Tests of Three-Dimensional Models in Calm Water and Waves,’ NSRDC Report 4095 (Sep 1973). 3. Martin, M., “Theoretical Predictions of Motions of High-Speed Planing Boats in Waves,’’ DINSRDC Report 76-0069 (Apr 1976). 4. Shuford, S.L., Jr., ““A Theoretical and Experimental Study of Planing Surfaces Including Effects of Cross Section and Plan Form,” National Advisory Committee for Aeronautics Report 1355 (1957). 5. Fridsma, G., “A Systematic Study of the Rough-Water Performance of Planing Boats,’’ Davidson Laboratory, Stevens Institute of Technology Report R1275 (Nov 1969). 33 hie bins vi bay eG i xP aU a ookrureD: APPENDIX A EVALUATION OF HYDRODYNAMIC FORCE AND MOMENT INTEGRALS The hydrodynamic force the craft experiences in the vertical direction as derived in the text is: ‘ dm,V ; ; F,=-f {mV-u 9E MV iste IDV }cosoat+ J apeadt where U = XGG cos @ - (z-w,) sin 0 and V = Xcg sin ®@ + (z-w,) cos 8 -6§ Another force acting in the vertical direction is the weight of the craft. The first two terms of the integral are evaluated by making the substitutions V =Xeq sin? —-@€ + Z—, cos'8-w, cos @ cE) ars cos 6 - ZeG sin 0) +w, 6 sin 0 . Ow a = -8 - = cos 8 A ONE) ; = aa o£ 0& dw ow aah z Taha he ouRE and noting that dUV = oe on Jom cine Using the previously described substitutions, the force becomes am, fov BE d§=-UVm, 35 cos 0) aoe {-M, cos 8 Z., - M, sin 6 nee + Q,6 + M, 0 (Zeg sin 6 - ets dw, A + fm, a cos 0 ag- ['mw,0 sin 6 dé Q Q dw, vi dw, = V — si + —_ fm aE sin 6 dé A m,U DE cos 6 dé -UV - | Vm,dé - 2 mM, fight i m, dé pf ca,cb¥ dé} cos 6 + [ apenct Q where M, = [mae Q and Q, = [meat Q This is essentially the form in which the integrals have been computed in the program. The rate of change of the sectional added mass in the third term of the integral expression is derived by relating it to the rate of change of depth of fluid penetration of the section. The added mass of a section is assumed to be equal to ma aka me pb? for which the time derivative is m, = k, tpbb where b is the instantaneous half-beam of the section, and k, is an added-mass coefficient, assumed to be constant. A value of k, = 1.0 was used in the computations contained in this report. For sections with constant deadrise, which is an imposed limitation of this work, the half-beam is related to the depth of penetration by b=dcotB 36 where d is depth of penetration, and @ is deadrise angle. Taking into account the effect of water pileup, the effective depth of penetration d, is, according to Wagner d, =7/2d e and b=d, cot 6 = m/2 d cot B where 77/2 is the factor by which the wedge immersion is increased by the pileup. Using this expression for the half-beam, the rate of change of sectional added mass becomes m, = kampb(n/2 cot B)d This expression is valid for penetration of the section up to the chine. When the immersion exceeds the chine, the sectional added mass is assumed to be constant, i.e., 2 may = koma plOme rh, = 0 where Dee is the half-beam at chine. The submergence of a section in terms of the motions is given by h=z-r where z = Zo, - & sin 0+¢€ cos 0 ia COS (kGee +£ cos @ +f sin 9) + wt} For wavelengths which are long in comparison to the draft and for small wave slopes, the immersion of a section measured perpendicular to the baseline is approximately Dial ~ cos @ -ysin#@ where v = wave slope 37 The rate change of submergence d is given by z-Tf + (z -1r) d(cos 8 - usin 0) eo" ee cos 6 - usin 0 (cos @ - v sin 6)2 ot Since immersion (z-r) is always small in the valid range of the previously described expression, the relationship can be further simplified to WG gee ~ cos 6 -vsin 8 and (z -1) m, ~ k,tpb(n/2 cot 8) ~—— The expansion of the integral expression for the hydrodynamic moment in pitch follows the procedure used for the vertical force. The results are summarized as follows Fo = -1,0 + Q, cos O Zog - Q,4Z¢g sin 8 - kg cos 8) dw, z ~ fm, cose =? kat + [m6 sin 0 w, dé Q ! Q + [ vingeat + [ocybveede Q Q +m,UVé| + [m,vuee stern Q uf dw, st; m,V —— sin @éd us dw, - / m,U = cos0@éd J m,U Gp cos 0 Edt + [aoe cos 6 €dé Q The only additional moments are the buoyancy moments. All other moments are considered to be zero for the specific problem considered in this report. 38 APPENDIX B COMPUTER PROGRAM DESCRIPTIONS OVERVIEW The equations of motions developed in the previous sections of this report have been solved by means of digital computer programs. Two major programs have been developed: the first (MAIN) solves the equations of motion using the Runge-Kutta-Merson integration algorithm and generates time histories that are stored on the system disk. The second (PLTHSP) generates California Computer Products Company (CALCOMP) pen plots from the disk files. All programs were designed to operate on the Control Data Corporation computer system, located at the David W. Taylor Naval Ship Research and Development Center in Carderock, Md. Descriptions of input data required to execute the programs, job control cards, and programs follow. Sufficient detail is presented for this appendix to serve as a manual for use and maintenance. JOB CONTROL CARDS FOR PROGRAM MAIN Job control cards for program MAIN which computes time histories of the motion variables, are described as follows. If CALCOMP plots are not desired, TAPES need not be cataloged. Job Control Language Card: Comment Job Card Charge Card Standard facility card Standard facility card REQUEST,TAPE9,*PF. REQUEST,TAPE2,*PF. REQUEST,TAPE4,*PF. ATTACH,BINAR,SEFZARNICKNEWB, ID=XXXX. ATTACH,NSRDC. LDSET(LIB=NSRDC). BINAR. REWIND,TAPE2. REWIND,TAPE4. COPY(TAPE2,OUTPUT) COPY(TAPE4,OUTPUT) Bg, Reserves space for CALCOMP plot data Print output file 1 request Print output file 2 request Attaches binary run file Attaches library routines Loads library routines Loads and executes run file Rewinds time-history files for printing Prints time-history file Prints time-history file Job Control Language Card: Comment CATALOG,TAPE9, SEFZARNICKDATA.., Catalogues file for plot. ID=XXXX. (SEFZARNICKDATA CAN BE ANY NAME) 7/8/9 END OF RECORD DATA CARDS (1-5) 6/7/8/9 END OF FILE INPUT DATA CARDS FOR PROGRAM MAIN Input data used by program MAIN are read from data cards in NAMELIST and in standard format. A description of the FORTRAN symbols appearing in NAMELIST follows. For simplicity in the text that follows, it is assumed that NAMELIST input occupies only one card. More cards can be used if necessary. Card 1(NAMELIST FORMAT, / i) A The absolute error for KUTMER (six values) NPRINT If=1, print normal output If=2, matrix, inverse matrix, F-column matrix, and KUTMER results If=3, integral results If=4, calculated values constant for given input values NPLOT If=0, no plot If=1, printer plot of results END Number of runs to be made W Weight of craft in pounds BL Boat length in feet TZ Thrust component in z direction TX Thrust component in x direction XECG Distance from center of gravity to center of pressure for drag force in feet XP Moment arm of propeller thrust XD Distance from center of gravity to center DRAG Friction for drag force RO Wave height LAMBDA Wavelength RG Radius of gyration in feet a Propeller thrust in pounds GAMMA Propeller thrust angle in degrees 40 Card 1 (continued) ECG Longitudinal center of gravity NCG Vertical center of gravity, nondimensionalized by ship length KAR Added-mass coefficient BETA() Dead-rise angle in degrees EST() Station position in feet NUM Number of stations XA Initial time XE Stop time HMIN Minimum step size HMAX Maximum step size EPS Error criterion Card 2 (Format 8F10.0) (X(1),I=1,6) Initial conditions X(1) Velocity X(2) Z X(3) i) X(4) xX X(5) Z X(6) 0 deyrees Card 3 (8F10.0) START Time to turn on (RMP) function (see page 48) RISE Duration of RMP Card 4 (8F10.0) TME Time at which integration interval is to be changed* HMX New maximum interval size after TME HMN New minimum interval size for KUTMER to subdivide *If this option is not used set TME to stop time on run. 4] Card 5 (8F10.0) PERCNT Percentage of boat length subtracted from longitudinal center of gravity to obtain X - point where acceleration computations are made JOB CONTROL CARDS FOR PROGRAM PLTHSP Job control cards for program PLTHSP which generates CALCOMP plots of time histories computed by program MAIN are described in this section. Job Control Language Card: Job Card Charge Card REQUEST,TAPE7,HI. Comment Standard facility card Standard facility card Tape for CALCOMP plot data VSN(TAPE7=CK0323). Volume serial number of tape for CALCOMP plot ATTACH,CALC936. Attaches CALCOMP library routine ATTACH, BINAR,SEFZARNICKPLOTB, Attaches plot program run file ID=XXXX. LDSET(LIB=CALC936) Loads CALCOMP library routines BINAR. Runs plot program 7/8/9 END OF RECORD DATA CARDS 6/7/8/9 END OF FILE INPUT DATA CARDS FOR PROGRAM PLTHSP Two or three data cards are made ready by PLTHSP, depending on the options selected. Standard input format is employed. A description of the necessary data cards follows. Card 1 (8F10.0 Format) XAXIS Length of x axis in inches YAXISP Height of pitch component axis in inches YAXISH Height of heave component axis in inches latit Height of lettering in inches Card 2 (110 Format) IA If=0, no plots for bow acceleration and center of gravity acceleration If=1, plots previously mentioned information 42 Card 3 (8F10.0 Format) - Only Necessary If |A=1. YAXISB Height of bow acceleration axis in inches YAXISC Height of CG acceleration axis in inches PROGRAM MAIN Program MAIN reads all necessary input data from cards, sets up initial values, computes constants, calls KUTMER to determine the state variables at TIME for the period from XA to XE in increments of HMAX. A table state variables is created for every PTIME-th value. The values for \/H and 0/27 H/d are calculated and printed. If the plot option is on, a printer plot will be produced. Subroutine COMPUT(X) This routine computes pitch moment NL and lift force FL, excluding added mass terms, using values of integrals computed in subroutine FUNCT. The argument X contains the state vector. Subroutine DAUX This subroutine is called from KUTMER or EULER. It determines the values of m,. b, and bl*, based on the following equations hed) = Zeg = E(I) sin 8 + (1) cos 6 - r(1) where r(I) = r, cos k ixee + €(1) cos 0 + €(I) sin @ + ct] Then for hy (1) 0) hy, OD) = ae at) Ga where V(I) = -rk sin 0 ixee + &(1) cos 6 + (1) sin 6 + ct] If d(1) > b,,() tan (BC) 2/n) 43 set m, (1) * Mamax 1) b(1) = b,, WD b1(1) =0 Mamax(l) = k(1) (p/2)mb? (1) If d(1) 4X97HCG ACCLo//) WRITE (4992) TIME» (X(I) sIT=1 96) sNL FL »BWACL »CGACL WRITE(9) TIME» (X(1) 5 1=496) »BWACL »CGACL KOUNT = KOUNT?1 FX(101B)=X(5) FX (291B) =X (6) IKUTM= (XE=xA) /HMAX*.05 IKUTM = (TME=XA) ZHMAX @ (AE=TME) SHMX @ .05 FIRST=0.0 NEQS=6 IKUTS=0 START OF INTEGRATION LOUP CONTINUE NPRINT = IPRINT # # @ # # ® CHECK PITCH .GT, 5236 RADIANS IF (X(6) eGT.05236)G0 TO 853 ® ® @ # © » PERFORM INTEGRATIONS IF (TIME LT. TMeOReTMEEQeXE) GU TO 98 IF (IPT.EQ.1) GO TO 98 HMIN = HMN HMAX = HMX FIRST = 0.0 CONTINUE CALL KUTME9 (NEQS 9 TIME »HMAX 9X 9EPSE sAgHMINgF IRST) IKUTS=IKUTS¢1 IF (FIRST.E9.2)GO0 TO 861 IF (KOUNT NF el eANDeKOUNTeNE 041) GO TO 99 WRITE (4991) KOUNT=1 *# # @ # # # WRITE OUT TIME INTERVAL RESULTS WRITE (4992) TIME» (X(I) 9 I=] 96) s9NLoF lL oBWACL »CGACL WRITE (6993) Tl oT2oT39T4eTSoT65T7 9 T898MMoBF WRITE(9) TIMEs (X(1) 9 1=496) »BWACL »CGACL IF (TIME sLT.TMceOReoTMEcEQeXE) GU TO 200 IF (IPT.EQe1) GO TO 200 CALL PLUTE? (FX5XAgHMAX »LAMBDA sIBoNWAVEs IPT) IPT = 1 52 MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 NAAANANANANANAANAANAANAANAINAND 20 CONTINUE 1821861 FX(1,1B)=x(5) FX(2s1B) =x (6) 93 FORMAT(" ",10E10.4) 92 FORMAT (1Xo11(F10.402Xx)) 190 CONTINUE ~ KOUNT=SKOUNT ¢1 IF (NWAVE.GT.0)GO TO 21 IF (TIME GT. WAVE) NWAVE=KOUNT 21 CONTINUE IF (TIME LE. XE AND. IKUTSoLT.IKUTM)GU TO 851 WRITE (29852) 854 CONTINUE 852 FORMAT( “ END OF KUTMER") 853 CONTINUE CALL PLUTE2 (FX »XAoHMAX »LAMBDA s IB yNWAVE 5 IPT) # # @ @ # @ # CHECK FOR LAST RUN IF NOT CYCLE HACK TO READ NEW DATA FOR NEXT RUN IF (ENO.NE.1)GU0 TO 8 GO TO 999 # # @ # KUTMER ERROR MESSAGES 861 WRITE (69862) 862 FORMAT ("' ERROR CRITERION IN KUTMER CAN NOT BE MET") WRITE (6956) (X(I) 9I=196) WRITE (6986) TIME 86 FORMAT (" TIME ='"9F10.4) IF (END.NE.1)GO TO 8 GO TO 853 999 CONTINUE END FILE 9 END SUBROUTINE PLUT2 (F oF MINoF MAX gNVAR »NFUN NI 9N9 XO 9DELX) PLUT FIRST N POINTS OF UP TO 26 FUNCTIONS F(X) F(I9J) CONTAINS THE VALUE FOR THE JTH POINT OF THE ITH FUNCTION FMIN(I) AND FMAX(I) CONTAIN THE MIN AND MAX ORDINATE VALUES FOR THE ITH FUNCTION, NVAR(I) AN ARRAY OF TITLES FUR THE VARIOUS FUNCTIONS TO BE PLOTTED AGAINST THE ABSCISSA NF UN NUMBER OF FUNCTIONS TU BE PLOTTED - OIMENSION OF NVAR» FMIN» FMAX Nl USED ONLY IN F(Nl»l) AS PASSED DIMENSION N NUMBER OF POINTS IN A SINGLE PLOT FRAME x0 FIRST ABSCISSA VALUE DELX ABSCISSA INCREMENT 1 2 ro DIMENSION 2STEP (26) oF (NI oN) oF MIN (NFUN) oF MAX (NFUN) » VLAST (26) » VF 19ST (26) HEAD (6) »STEP (26) INTEGER CH(26) »NVAR( NFUN) »DOT»ASTERsPLUS » BLANK INTEGER C INTEGER A (101) DATA BLANK ,DOT,ASTER»PLUS/IH olHeolH®slH+e/ DATA CH(1) »CH(2) »CH(3) oCH (4) »CH(S) »CH(6) »CH(7) »9CH(8) »CH(9) »CH(10) 7 HA 4» LHB » LHC » LHD » LHE » LHF » 14G » IHH oIHI olHJ # DATA CH(11) »CH(12) »CH(13) 9CH(14) »CH(15) »CH(16) »CH(17) sCH(18) f 7 AK 4 LAL » JAM » JAN » JHO » LHP » ILH@ » ILHRS DATA CH(19) »CH(20) oCH(21) »CH (22) »CH(23) »CH (24) »CH(25) »CH (26) 53 MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN MAIN PLOT2 PLOT2 PLOT2 PLoT2 PLOT2 PLOT2 PLOT2 PLOT2 PLOT2 PLOoT2 PLoT2 PLOT2 PLOT2 PLoT2 PLoT2 PLOT2 PLOT2 PLOT2 PLOT2 PLoT2 PLOT2 PLOT2 PLOT2 PLOT2 PLoT2 PLOT2 PLoT2 2 4 HS » HT » IHU » LHV 5» LHW 9» IF (NFUN.LE.0.OR.N.LE.0) RETURN C PRINT HEADINGS, 46 30 WRITE (6946) FORMAT (///) DO 40 I=l 5NFUN TENM=A0S (FMAX(1) 5ND KUTMER65 Y2(1) = YO(I) *HC/6.#F0(1) + (20/36) *HC#F2 (1) + (HC/6.) #F 1 (T) KUTMER66 IF (NPRINT .£Q.5) WRITE (69400) Y20T KUTMER67 INC = 0 KUTMER68 ese eer eee CHECK ERROR CRITERIA KUTMER69 DO 110 1=1.NO KUTMER7O ZZZ = ABS(YIL(I))=-A(I) KUTMERT1 IF (222) 85587,87 KUTMERT2 soso = 2 = = ABSOLUTE ERRUR KUTMER73 ERROR = ABS(.2#(Y1(I)-Y2(1))) KUTMERT4 IF (ERROR=A(1)) 100,100,590 KUTMERTS cee cee eee RELATIVE ERKOR KUTMER76 ERROR = ABS(e2-e2"Y2(1)/VI(1)) KUTMERT7 IF (ERROR=EPSE(I)) 1009100990 KUTMER78 Fe a ie Sule SINCE ERROR eGT. ERROR CRITERIA CHECK IF HC.GT H/KUTMERT9 ee eee ee IF YES THEN HALVE INTERVAL. OTHERWISE STOP, KUTMERBO X = 128.#AS (HC) =ABS (H) KUTMERB1 IF (X) 91,95,95 KUTMERB2 -2f2ee-+-ee ERROR TOU LARGE KUTMER83 WRITE (6992) I9T sERROR®HC KUTMER84 FORMAT (/18H FOR EQUATION NO. I2927Hs THE RELATIVE ERROR AT T = KUTMERBS e E15e8>5 4H IS 9£15.8913H STEP SIZE = 5615-8) KUTMER86 FIRST = 2. KUTMERB7 RETURN KUTMER88 =PHALF DAUX 19 COMMON /IN/ BM (120) ,81 (120) »VELIN DAUX 20 COMMON/UUT /NPRINT »NPLOT sEND VAUX 21 COMMON /SEAWAVE/ START »RISE ,RAMP DAUX 22 COMMON /WAVE/ RoPT(120) oZMAgZWMA EMAS 5 ZZWMA 9 ZWEMA 9 ZOWMA 9ECMAZ 9 DAUX 23 6 ZwOOT (120) DAUX 24 Cc DAUX 25 RAMP = RMP(TIME START RISE) DAUX 26 PIH = PI/2, DAUX 27 CT = C#TIME DAUX 28 Cx6 = CUS(x(6)) DAUX 29 SK6 = SIN(X(6)) DAUX 30 CaeeeeeeSET VALUFS UF MA AND B DAUX 31 00 75 I=1,NUM DAUX 32 PT(L) = (X(4) E (1) #CX6eN(1) #SA66CT) HK DAUX 33 R(I) = RU®COS(PT(I)) #RAMP DAUX 34 Cc # * # # # # # # COMPUTE HW SUBMERGENCE OF A POINT AND R- THE WAVE DAUX 35 Cc HW(T) IS IN THE FIXED COORUINATE SYSTEM DAUX 36 Si Cc 65 Cc Cc Cc Cc Cc 70 75 74 76 77 78 79 80 si 82 85 Cc Cc # # Cc 15 17 Cc # * 18 C # # HW(I) = X(5)-E (1) #SX6eN(1) #CX6-R (1) IF (HW(1) GT.) GO TO 65 CRAFT IS NOT SUBMERGED MA(I) = 0. Bl(1I)=0. B(I) = 0. GO TO 75 V(I) -RO#K*SIN(PT (1) ) *RAMP 0 (1) Hw(I)7(CA6—-V (1) #SX6) D(I) 1S IN THE BODY AAIS SYSTEM AND IS THE SUBMERGENCE IF (D(1) GE.TEST(I)) GO TU 70 CRAFT IS PARTLY SUBMERGED B(I) = O(1)#(1./TA)#PIH BL(T> D(T)#(1-/TA) #PI MA(I) KAR#PHALF #B (I) #B(1) GO TO 75 CHINE IS IMMERSED Bl ARRAY IS USED FOR THE INTEGRALS OVER THE PORTION OF THE HULL FOR WHICH THE CHINE IS NOT IMMERSED MA(I)=MMAX(1) B(I)=BM(1) B1l(I)=06 CONT INU IF (NPRINT.LTe%) GO TO 85 WRITE (6974) TIME FORMAT(" TIME = "yF10.4) WRITE (6976) (X(T) » T=1 96) WRITE (6977) (R(I) » I=] NUM) WRITE (6978) (HW(I) »I=1 NUM) WRITE (6079) ( B(I) »I=l NUM) WRITE (658%) (V(1I) o I= »9NUM) WRITE (6581) (O(1) » l=] NUM) WRITE (6982) (MA(I) » T=) »NUM) FORMAT(" X(I) "“96(2X0E120c6)) FORMAT (* R(T) %910F 10.4) FORMAT ( HW(I)%510F10.4) FORMAT (" 8(1)%510F10.4) FORMAT (* v(I)%510F10.4) FORMAT (" (I) %910F 10.4) FORMAT (" MA(I) “slOF1064) CONTINUE * @ # # # » COMPUTES NL AND FL ANU THE ASSOCIATED INTERGALS CALL FUNCT (X) IF (NPRINT.LT2%)GO TO 17 WRITE (6915) TXoFLoDRAGoTZoWgNboXDoT oXP FORMAT(" ",10E12.6) CONTINUE # # # # # COMPUTE THE F VECTOR F (lol) = \TX#FL#SX6=DRAG#CX6 F(151)=0.0 F(20l) = T2*FL#CX6*DRAG*SX6eW F (391) =NL=DRAG#XD+T#xXP IF (NPRINT.LT.3)G0 TO 18 WRITE (6910) (F (Tol) »I=193) CONTINUE # &# # # # COMPUTE THE A MATRIX A(lol) = M+MASS#SX6#SX6 58 DAUX DAUX DAUX DAUX DAUX DAUX DAUXK DAUX DAUX DAUX DAUX DAUX VAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX OAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX A(lo2) = MASS*®SKX6#CX6 A(leo3) = =)A*#SX6 A(ls2) = oe A(l9o3) = Al2el)= AIRE) A(202) = M+MASS#CX6#CX6 A(293) = -QA#CX6 A(391)=A(1,.3) A(392)=A(2.3) A(393)=SIT*IA IF (NPRINT.LT2£3)GO TO 25 WRITE (6012) (A(I 91) oI=193) WRITE (6913) (A(I 92) oL=1 93) WRITE (6914) (A(1,3) »I=1,3) # @ # @ # INVERT THE A MATRIX CALL MATINS (A93939F 91919DETERMsIDo INDEX) IF (ID.EQ.2) WRITE (6026) FORMAT (*! MATRIX IS SINGULAR my) Ceeee#eA ON RETURN WILL CONTAIN THE INVERSE MATRIX c*#* 25 26 Cc Cc Cc c## 10 12 13 14% 39 35 40 Cc I0=2 MATRIX IS SINGULAR =1 INVERSE WAS FOUND ® # # © # COMPUTE THE RIGHT HAND SIDE RHS(1) = F(lsl) RHS(2) = F(291) RHS(3) = F(391) RHS(1) = 0.0 RHS(4) = X(1) RHS(5) = X(2) RHS(6) = X(3) FORMAT (" F(T oh) %93(2X5E12.4)) FORMAT ("' A(Tol) "s3(2XsE12.4)) FORMAT (" A(T92) "s3(2XsE12.4%)) FORMAT ("" A(1593) "93(2XsE1264)) IF (NPRINT.LT.2) GO TO 40 WRITE (6912) (A(Io1) ,IT=193) WRITE (6013) (A(1,92) »I=1 93) WRITE (6914) (A(I93) 9I=193) WRITE (6935) (RHS(I1) »I=1 96) FORMAT ("* RHS(I) %96(2X9E1206)) CONTINUE RETURN END SUBROUTINE FUNCT (X) REAL KAR REAL IAs TAAgIPART 9K oKPI o>MAgMASS NL »NCGo IT oMoMMAX oN INTEGER END DIMENSION IPART (120) C1 (120) »C2(120) 5 26 (120) oZ7(120) 0X (6) 9 VMAA (120) e@ee@fFf 11 (120) 902 (120) 903 (120) 904(120) 505 (120) 906(120) » QPART (120) 9Z1 (120) 922 (126) 923 (120) 924 (120) 925 (120) » DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX DAUX FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT COMMON /SHIP/ MASSsCINT sQA9CEsCE22CE3sDMU sEDMUsE2DMU 9 E 3DMUs BF s BMM oF UNCT NL oFLoTAsE (120) * COMMON JCONST/ NCGoECGoPI »DPR»oRPD »sGRAVTY »RHO 9K »NUM9MA (120) »CDo TAs e B(120) sBETAsHW (120) 9» TZ9DRAGgWoXDoTs XPoMolTo DELTAS o TXoEST (120) sCoROsKARgMMAX(1 0) op TEST (120) » e N(120) »>PHALF 59 FUNCT FUNCT FUNCT FUNCT FUNCT & 50 si 60 COMMON /INZ BM(120) oB1 (120) »VELIN COMMON /U0UT /NPRINT »NPLOT »END FUNCT FUNCT COMMON /WAVE/ R(120) sPT (120) o> 2MAgZWMAgEMAS,ZZWMA, ZWEMA, Z2WMA,E2MAZF UNCT 9ZWOOT (120) e COMMON /INTER/ II oKTT (10) oDIFF (10) COMMON /SEAWAVE/ START »RISE,RAMP COMMON /TEST/ VMA #2 2 #@ #@ » # @ INITIALIZE INTEGRAL SUMS 0.0 (1) ®SIN (XK (6)) oX (2) @CUS (X (6) ) Sx6 = SIN(x(6)) CX6 = COS(X(6)) WO = K#C # @ # @ # # SET UP THE FUNCTIUNS FOR THE INTEGRALS DO 90 I=1l»NUM IPART (1) =E (1) #E (1) #MA(T) QPART (I) =E (1) #MA(T) ZWOOT(I) = -RU#WO#SIN(PT (I) ) #RAMP U = X (1) #CX6—-X (2) *5X6e¢ZWDOT (1) #SX6 VEL = VPART=X(3) #E (1) -ZwDOT (1) #CX6 N N a = > 1 0 tb Z1(1T) = MA(I)*#ZWOOT(T) Z2(1) = -Ma(1I)*#COS(PT(I)) #RAMP Z3(1) = E(1) #221) 24(1) = E(T) #2101) Z5(1) = U#z2(1) Z6(1) = E(7)#25(1) Z7(1) = MA(I) #VEL#U IF (VELeLE.0.) GO TO 60 IF (BI(1) LE. 0.0) GO TO 50 DRDT = ZwOOT (1) #(X(L) Cox (3) ®(N( 1) @CX6=E (1) #SX6) ) /C DI(1) = — VEL®BL (1) #(X (2) x (3) # (CXO*E (1) #SX6#N(T)) GO TO 51 DI(1) = 0. CONTINUE p2¢1) = €(1)*D1(1) Cl(1) = VEL#VEL#B(I) C2(I) = €(1)#C1 (1) GO TO 6 OAS) en Oa D2(I) = 0. C1(1) = 06 C2(1) = 0. 60 FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT (PAGE 4 OF NOFUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT =-DRDT) FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT FUNCT C Cc 6l CONTINUE 03(1) 04 (1) PIH = -OS(I) 06(T) = Z22(1) *VEL = E(1)#03(T) PIv2, = B(I)*(HwW(1)-B(I) #TAs20) = 05(1) *E(1)*.5 CONTINUE RHOG=RHU#GRAVTY # #® # # ® SET UP THE FUNCTIONS FUR THE INTEGRALS (PAGE 5 UF NOTES)FUNCT 85 PIH = KPI = PI7/2, KAR#PI EVALUATE INTEGRALS USING TRAP METHOD 91 93 9% I = 1 INDEX CALL CALL CALL CALL CALL CALL CALL CALL CALL CALL CALL CALL CALL CALL CALL CALL CALL CALL CONTI MASS QA = OMU EDMU E2DMu E30MuU BF = BMM = ZMA = ZwWMA EMAS ZZWMA ZWEMA Z2WMA E2MAZ CONTI IF ( INDEX I2zI GO TO = } TRAP (MACINDEX) »DIFF (1) »KTT (1) » TMASS) TRAP (QPART (INDEX) sOIFF (I) sKTT (1) 9QAl) TRAP (C1 (INDEX) oOIFF (1) oKTT (1) CEA) TRAP (C2 (INDEX) »DIFF (1) oKTT (1) »CE2A) TRAP (IPART (INDEX) sDIFF (1) sKTT(1I) 9 IAA) TRAP (D1 (INDEX) »OIFF (1) 9K TT (1) 9OMUA) TRAP (D2 (INDEX) »DOIFF (1) »KTT (1) 9EDMUA) TRAP (N3 (INDEX) »DIFF (1) »KTT( 1) 9E2DMUA) TRAP (N& (INDEX) oDIFF (1) »9KTT (1) s9E3DMUA) TRAP (DS (INDEX) »DIFF (1) oKTT(I) 9 BFA) TRAP (N6 (INDEX) oDIFF (1) »KTT (1) »BMMA) TRAP (Z1 (INDEX) oDIFF (1) sKTT (1) 9 ZMAA) TRAP (72 (INDEX) oDIFF (1) oKTT (1) 9 ZWMAA) TRAP (73 (INDEX) »OIFF (1) oKTT(1) 9EMASA) TRAP (74( INDEX) »DIFF (1) »KTT(1) 9 ZZWMAA) TRAP (75 (INDEX) »DIFF (1) »9KTT (1) 9 ZWEMAA) TRAP (76( INDEX) »DIFF (1) »KTT(1) »9Z2WwMAA) TRAP (77 (INDEX) »DIFF (I) »KTT (I) s9ECMAZA) NUE = MASS ¢« TMASS QA + QAl IA + TAA CE + CEA cE2 + CE2A OMU + DMUA = EDM! + EDMUA = E20MU «+ E20MUA = E30MU «+ E3DMUA BF + CHOG*BFA BMM + RHOG#®BMMA ZMA+ZMAA = ZwMA+ZWMAA = EMAS+EMASA ZZwMA+ZZWMAA ZWEMA+Z2WEMAA Z2wMA+Z2WMAA E2MAZ+E2MAZA NUE 1.GE.11)GO 10 92 2 INDEX+KTT (I) =} ol 91 92 CONTINUE 61 FUNCT 77 FUNCT 78 FUNCT 79 FUNCT 80 FUNCT 81 FUNCT 82 FUNCT 83 FUNCT 84 FUNCT 86 FUNCT 87 FUNCT 88 FUNCT 89 FUNCT 90 FUNCT 91 FUNCT 92 FUNCT 93 FUNCT 94 FUNCT 95 FUNCT 96 FUNCT 97 FUNCT 98 FUNCT 99 FUNCT100 FUNCT101 FUNCT102 FUNCT103 FUNCT104 FUNCT105 FUNCT106 FUNCT107 FUNCT108 FUNCT109 FUNCT110 FUNCT111 FUNCT11l2 FUNCT113 FUNCT114 FUNCT115 FUNCT116 FUNCT117 FUNCT118 FUNCT119 FUNCT120 FUNCT121 FUNCT122 FUNCT123 FUNCT124 FUNCT12S FUNCT126 FUNCT127 FUNCT128 FUNCT129 FUNCT130 FUNCT131 FUNCT132 FUNCT133 FUNCT134 FUNCT135 Cc # * @ #@ # # # CALL COMPUT TO FIND THE VALUE OF NL AND FL USING FUNCT136 Cc THE VALUES OF THE ABOVE INTEGRALS FUNCT137 CALL COMPUT (xX) FUNCT138 Cc FUNCT139 IF (NPRINT.LT.3) GO TO 11} FUNCT140 IF (NPRINT.FEQe3) GO TO 108 FUNCT141 IF (NPRINT.EQ.4)GO TO 108 FUNCT142 WRITE (6997) (IPART(I) 9121 oNUM) FUNCT143 WRITE (6998) (QPART(I) »I21 9NUM) FUNCT144 WRITE(6999) (C1(I) 9 IT=1 9NUM) FUNCT145 WRITE (69100) (C2(1I) »I=l NUM) FUNCT146 WRITE (69161) (C3(I) »IT=2 oNUM) FUNCT147 WRITE (69102) (Ol (1) 9 I=] 9NUM) FUNCT148 WRITE (69103) (D2 (1) »J=1 NUM) FUNCT149 WRITE (69104) (03(1) 9 I=19NUM) FUNCT150 WRITE(6910S5) (D04(1) o2=) NUM) FUNCT151 WRITE (69106) (05(1) 5 I=] NUM) FUNCT152 WRITE(69112) (D6(1) 9 I=1 9NUM) FUNCT153 WRI FE (69113) (21 (1) » T=1 9NUM) FUNCT154 WRITE (69114) (Z2(1) 9 T=] 9NUM) FUNCT1S55 WRITE (60115) (23 (1) 9 L=1 »NUM) FUNCT156 WRITE (69116) (24(1) 9 T=1 9NUM) FUNCT157 WRI FE (69118) (25 (1) 9 T= oNUM) FUNCT158 WRITE (69119) (26 (1) 9 T=] 9NUM) FUNCT1S9 WRITE (69120) (27 (1) » T= 9NUM) FUNCT160 WRITE (69107) KPI »RHOG,PIH FUNCT161 108 WRITE (69109) MASSsCINT»QA9CE,CE2,CE3 FUNCT162 WRITE (69121)1A FUNCT163 121 FORMAT(* IA #9E10.4) FUNCT164 WRITE (69110) OMUsEDMU» E20MU 9 E 3UMU 9 BF 9 BMM FUNCT165 WRITE (69117) ZMAyZWMA sEMAS 9 ZZWMA, ZWEMA 9 ZZWMA gE2MAZ FUNCT166 Cc ## # @ &@ & @ © @ © FORMATS # #@ #e #@ toe ee ee FUNCT167 96 FORMAT (" CPART(1)%510(2XE104)) FUNCT168 97 FORMAT(" IPART (I) %910(2x9E1004)) FUNCT169 98 FORMAT(" OPART (I) %510(2X0E100¢4)) FUNCT170 99 FORMAT(" Cl %910(2X%9E1004)) FUNCT171 1@0 FORMAT(" C2 %910(2X9E10.4)) FUNCT172 101 FORMAT(" C3 %510(2XsE10.4)) FUNCT173 102 FORMAT(" Ol %510(2X2E1004)) FUNCT174 103 FORMAT(" p02 %o10(2X9E10¢4) ) FUNCT175 104 FORMAT(" D3 5 10(2X9E1004) ) FUNCT176 105 FORMAT(" 04 %910(2X%9E1004)) FUNCT177 106 FORMAT(" DS %510(2X9E1004)) FUNCT178 122 FORMAT(" 06 %910(2x%9E1004)) FUNCT179 107 FORMAT(" KPHI "sE100%9"RHUG "E1004" PHIH 961004) FUNCT180 109 FORMAT(" MASS "“sE10c49" CINT "“sEL0c4%9" QA “sE10c49" CE "sE100%s FUNCT181 ~ @0CE2 5610.40" CE3 9E10.4) FUNCT182 110 FORMAT (" DMU "9610.45 EDMU % 9610049" E2DMU 9610049" E3DMU ty FUNCT183 #E10.45" BF 5E10.45" BMM "5E10.4) FUNCT184 113 FORMAT(SH 21 910(2X5E10.4)) FUNCT185 114 FORMAT(GH 22 910(2X%sE1004)) FUNCT186 115 FORMAT(4H Z3 910(2X%sE1004)) FUNCT187 116 FORMAT(4H 24 910(2X5E10-4)) FUNCT188 118 FORMAT(4H 25 910(2X5,E10.4)) FUNCT189 119 FORMAT(4H 26 o10(2XxX5E10.4)) FUNCT190 120 FORMAT (4H 27 910(2X5E10.4)) FUNCT191 117 FORMAT(SH 2MA 9€10.45,6H ZWMA 9610.04 96H EMAS 961004 FUNCT192 es TH Z2WMA 5E10.4 97H ZWEMA »5E10.497H ZZWMA 9E10.4, FUNCT193 S 7H E2MAZ 5£10.4) FUNCT194 62 lil 20 CONTINUE RETURN END SUBROUTINE COMPUT (xX) OIMENSION x (6) REAL KAR KPI REAL NLoMASSoNCGoMoIToITAagK oMAoMMAX oN INTEGER END COMMON /SHIP/ MASS sCINTsGAsCE sCE2CE3, DMU »EDMU,E2DMU 9 E 3DMU 9 BF » BMM, COMPUT NL oFLoIAsE (120) ~ COMMON JCONST/ NCGeECGoPI sDPReRPD »GRAVTY oRHO 9K »NUM MA (120) oCDo TAs ry B(120) sBETAgHW(120) 5 TZ 9>DRAGoWoXDoToXPoMolTo ¢ DELTAS + TXoEST (120) 6CoROsKARyMMAX (1 0) sTEST(120) » N (120) »PHALF *COMMON/OUT /NPRINTs NPLOT »END COMMON /TEPMS/ TloT20T30T4%oTS9T69T79T8 COMMON /WAVE/ R(120) »PT (120) » ZMAy ZWMAgEMAS 5 ZZWMA 9 ZWEMA y ZOWMA g E2MAZ»ZWDOT (120) * COMMON JTEST/ VMA Cx6 = COUS(x(6)) SxX6 = SIN(X(6)) wO = KeC CONS] = RO#WO#WOFCX6 CONS2 = (KPT#RHO#PTH/TA) /CX6 CONS3 = RO#WOFKECX6eSKE CONSS = RO#WOFK#®CX6#CX6 TERMI = X(1) #®CX6 TERM2 = X(2)#SX6 UVNUM = (X(1) #CX6=(X (2) —ZWDOT (NUM) ) #SX6) # © (X (1) ®SX6—X (3) #E (NUM) + (X (2) -ZWDOT (NUM) ) #CX6) ZMA = ZMA#K (3) #SK6 ZZWMA = ZZWMA*X (3) #SX6 ZWMA = ZwWMA#®CUNS1 EMAS = EMAS#CONS] DMU = DMU*CONS2 EDMU = EDMU#CONS2 CE = CE*CD#RHO CE2 = CE2#cD#RHO E2DMU = E2NMU*CONS3 E30MU = E3DMU*CONS3 ZWEMA = ZWEMA®CUNS4 Z2WMA = Z2WMA*CONS4 Tl = QA#X (3) #(TERMI-TERM2) Tl = Tl + 2ZWMA = EMAS T2 = EDMU T3 = CE2 T4 = MA(NUM) ®E (NUM) ®UVNUM + E2MAZ + E3DMU - Z2WMA + BMM NL = Tl + T2 + T3 + TG + BMM TS = MASS®X (3) # (TERM2=TERM1) T5 = 15 + 2WMA = ZMA T6 = -DMU 17 = -CE 63 FUNCT19S FUNCT196 FUNCT197 COMPUT COMPUT COMPUT COMPUT COMPUT COMPUT oorouft wh COMPUT COMPUT10 ComMPUTI1 COMPUT 12 COMPUT13 COMPUT 14 COMPUT15S COMPUT16 COMPUT17 COMPUT18 COMPUT19 COMPUT20 ComPuT21 COMPUT22 COMPUT23 COMPUT 24 COMPUT2S COMPUT26 COMPUT27 COMPUT28 COMPUT29 COMPUT30 COMPUT31 COMPUT32 COMPUT33 COMPUT 34 COMPUT3S COMPUT 36 COMPUT37 COMPUT38 COMPUT39 COMPUT4O COMPUT41 COMPUT42 COMPUT43 COMPUT4&4 COMPUT4&S COMPUT46 COMPUT47 COMPUT48 COMPUT49 COMPUTSO COMPUTS1 COMPUTS2 COMPUTS3 COMPUTS4S COMPUTSS COMPUTS6 COMPUTS7 T8 = =-MA(NUM) #UVNUM = E2DMU + ZWEMA BF = BF/CX4 FL2=TS+¢T6+T7+T8-BF IF (NPRINT.LT.3)GO TO 30 25 CONTINUE WRITE (6910) NLoFL 10 FORMAT(" NL = ",E12.60" FL = 9612.6) 30 RETURN END SUBROUTINE INPUT C# # © @ © @ © DEFINITION OF INPUT VARIABLES XA = INITIAL TIME KE = FINAL TIME HMIN = MINIMUM STEP SIZE HMAX = MAXIMUM STEP SIZE EPSE = RELATIVE ERROR CRITERIUN USED FOR VALUES OF Y GT A EPs 2 ERROR CRITERION IN KUTMER A = ABSOLUTE ERROR CRITERIA USED IN KUTMER NPRINT = 1 FINAL PRINTOUT = 2 MATRIX INVERSE MATRIXoF COLUMN MATRIX»sAND KUTMER RESULTS 3 INTEGRAL VALUES 4 CALCULATED VALUES-CONSTANT FOR GIVEN INPUT VALUES —S it il NPLOT = NO PLOT = PRINTER PLOT END = NUMBER OF RUNS M = MASS OF CRAFT WwW 2 WEIGHT OF CRAFT Z = THRUST COMPONENT IN Z UIRECTION = THRUST COMPONENT IN x DIRECTION XECG = DISTANCE FROM CG TO CENTER OF PRESSURE FOR NORMAL FORCE xP = MOMENT ARM OF PROPELLER THRUST xD = DISTANCE FROM CG TO CENTER UF PRESSURE FOR DRAG FURCE KA(I)= ADDED MASS COEFFICIENT AN ARRAY GIVEN THE VALUE KAR wHICH IS READ IN BM(I) = BEAM AT FREE SURFACE UR AT CHINE DRAG = FRICTION DRAG K = WAVE NUMBER RO = WAVE HEIGHT NU = WAVE SLOPE NUM = NUMBER OF STATIONS BL = BOAT LENGTH LAMBDA = WAVE LENGTH RG = RADIUS OF GENERATION IN FEET T = PROPELLED THRUST IN LBS GAMMA = PROPELLER THRUST ANGLE IN VEGREES DELTAS=STATION SPACING IN FEET ECG = LONGITUDINAL CENTER OF GRAVITY NCG = VE?TICAL CG BETA(I) = DEAD RISE NO(I) = HEIGHT OF MEAN BUTTOCK RHO = DENSITY OF WATER GRAVTY = GoAVITY FT/SEC##2 DPR = DEGPEES PER RADIAN RPD = RADIANS PER DEGREE 3.14159 2 oo 0 o © DADMDAANAANAAANANANAANAMDDMAMNMNAMNANNANANNANANANANNANNANANAANAANANAANAANANANA 64 COMPUT58 COMPUTS9 COMPUT60 COMPUT61 COMPUT62 COMPUT63 COMPUT64 COMPUT65 COMPUT66 COMPUT67 COMPUT68 INPUT 2 INPUT 3 INPUT 4 INPUT 5 INPUT 6 INPUT 7 INPUT 8 INPUT 9 INPUT 10 INPUT 11 INPUT 12 INPUT 13 INPUT 14 INPUT 15 INPUT 16 INPUT 17 INPUT 18 INPUT 19 INPUT 20 INPUT 21 INPUT 22 INPUT 23 INPUT 24 INPUT 25 INPUT 26 INPUT 27 INPUT 28 INPUT 29 INPUT 30 INPUT 31 INPUT 32 INPUT 33 INPUT 34 INPUT 35 INPUT 36 INPUT 37 INPUT 38 INPUT 39 INPUT 40 INPUT 41 INPUT 42 INPUT 43 INPUT 44 INPUT 45 INPUT 46 INPUT 47 INPUT 48 INPUT 49 ANNNANAANANAANN a0 aa0 EST(I) = STATION POSITION START = RISE = eo * @ #@ © # # » IC OPTIONS REAL IT 9K LAMBDA. MoMA o>MMAX gNUsNoNCGgNO MASS oNL og LAgKAR COMMON /CONST/ NCGeECGoPI sDPRaRPD sGRAVTY »sRHO 9K »NUM9MA(120) »CDoTAg B(120) sBETAsHW (120) > TZyDRAGgWoXDoToXPoMolTo DELTAS 9 TX9 EST (120) »CoROoKAR 9 MMAX (170) » TEST (120) 9 & 6 C COMMON /IN/ BM(120) 981 (120) »VELIN COMMON /IN2/ NO(120) »XA9XE9HMAXsHMINSA (6) o0EPSE (6) sLAMBDA C) eeamgeenonedee IC(1l) =1 USE WAVE Z DISTANCE IN COMPUTING LIFT COMPONENT OF NL AND FL INTEGER END N(120) »PHALF * COMMON JSHIP/ MASS oCINT 9QAoCE 9 CE2 0CE3 yDMU»EDMU 9 E2DMU 9 E3DMU 9 BF 9 BMM y NL oFL eo ITAgE (120) COMMON/OUT /NPRINT »NPLOT 9 END START TIME OF THE RAMP FUNCTION FOR SEA WAVE DURATION OF THE RISE FROM ZERO TO ONE OF THE RAMP COMMON /ACCEL/ XACCL »BWACL »CGACL BL NAMELIST/HSP/A gNPRINT »NPLUT gEND pWo Kh 9 TZ9 TX oXECGoXP oXDy DRAG »RGo Ts GAMMA sECU »NCG 9 KAR» ROgLAMBDA pNUM,BETASEST 9X Ao XE »9HMINgHMAX EPS o VELIN DATA A bie 015200015.00001501,,.V001».00001/7 DATA NPRINT»NPLOTsEND/1 9191/7 DATA WoBLoTZoTX9XECGyXP9XD 9DRAG»RU»sLAMBDA 9RGo T 9 GAMMA 9 ECGeNCGoKAR /166 930759640009 00416922659 095629280 00% 2232590.091,0/ “DATA NUM »BETASEST 477520.05 0.00905 .031255 6062509 0093755 2125005 0156255 6187509621875, 025000 5.281259 631250 9 6343759 6375005 0406255 6437509 0468755 0500995531255 656250 9 0593759 6625005 6656255 2687596718759 0750005. 781255 681250 5 0843759687500» 690625 937509 69687591.000 1.06250 .1.1250051.1875091 62500091 .312551.3750091.4375, 1.50051 ,5625 91 662591 668759167591) 0812551087591 09375 52005 200625 5 2.125 926187592 025920312992 037592643715 9265924562592 06255 2 68759267592 081259208150 9209315930059 30062553012593,1875, 30250053.3125 9363759304375 030593056255 3.6255 30687543615 / DATA XAsXE ,HMINsHMAXsEPS /0,0920.090025y0l0015/ DATA VELIN /19,62/ # @# 9 # # # READ IN AND WRITE OUT KUTMER PARAMETERS AND PROGRAM OPTIONS READ (S,HSP) WRITE (69HSP) OO 10 I=194 EPSE(I) = EPS # # # # # # SET UP CONSTANTS PI = 3.141592653589 GRAVTY=32.18 DPR=57 29577951308 RPD=.01 7453292519 65 INPUT 50 INPUT Sl INPUT S2 INPUT 53 INPUT 54% INPUT 55 INPUT 56 INPUT 57 INPUT 58 INPUT 59 INPUT 60 INPUT 61 INPUT 62 INPUT 63 INPUT 64 INPUT 65 INPUT 66 INPUT 67 INPUT 68 INPUT 69 INPUT 70 INPUT 71 INPUT 72 INPUT 73 INPUT 74 INPUT 75 INPUT 76 INPUT 77 INPUT 78 INPUT 79 INPUT 80 INPUT 81 INPUT 82 INPUT 83 INPUT 84 INPUT 85 INPUT 86 INPUT 87 INPUT 88 INPUT 89 INPUT 90 INPUT 91 INPUT 92 INPUT 93 INPUT 94 INPUT 95 INPUT 96 INPUT 97 INPUT 98 INPUT 99 INPUT100 INPUT1O1L INPUT1O02 INPUT103 INPUT104 INPUT10S INPUT106 INPUT107 INPUT108 aa0 IF (EST (NUM) .LT23075) STOP 3 COMPUTE NO AND BM ARRAYS DO 32 I=1,NUM IF(EST(I) .GE.0.75) GO TO 30 NO (T) =-0046875# (100-SQRT (EST (1) /00375-(EST (1) 0075) ##2,0) ) BM (1)=.375#SQRT (1.0-(EST(1) 7. 75-10) ##200) GO TO 32 30 NO(I)=0.0 BM(I) = 0.375 32 CONTINUE Ceeeteee2COMPUTE CONSTANTS AND INITIALIZE ARRAYS ANANNQANAANANANANANAN oO Cc IN M=W/GRAVTY RHO=1.99 IT=M#RG*RG K = 2.%PT/LAMBDA C=SQRT (GRAVTY/K) NU=RO#K PHALF = (PI/2.)#RHO BETA = BETA#RPD CO = COS(BFTA) TA = TAN(BETA) DO 60 I=1.NUM E(I) = ECG-EST(I) N(I) = NCG*NO(I) MMAX(1) = KAR®PHALF®BM(I) #BM(1) TEST(I) = (20%BM(1)#TA)/PI1 CONTINUE END=END*1 RETURN END SUBROUTINE PLUTER(FX»XA»HMAX »LAMBUA, IB yNWAVE 9 IPT) PUT: FX A TwO DIMENSIONAL ARRAY CONTAINING PITCH AND HEAVE VALUES AT EACH TIME STEP KA INITIAL TIME HMA x TIME INTERVAL» PTIME*HMAX = INTERVAL BETWEEN FX VALUES LAMBDA WAVELENGTH USED IN CALCULATING PITCH AND HEAVE RATIOES IB NUMBER OF FX VALUES NWAVE START OF VALUES AFTER WAVE IS COMPLETELY ON REAL ITsK »LAMBDA 9MoMA oMMAX gNoNCG INTEGER END DIMENSION FX (29400) »FMIN(2) »FMAX (2) sNVAR (2) COMMON /CONST/ NCGsECGoPI sDPR»RPD»GRAVTY »RHO 9K »NUMo9MA (120) sCDo TAs e B(120) sBETAsHW(120) »TZsDRAGoWoXDoToXPoMoITs DELTAS 9 TX oEST (120) sCoROoKAgMMAX (120) sTEST (120) 5 N(120) »>PHALF @ COMMON/OUT /NPRINT sNPLOT »END C# # # # ¢ # # # SET UP VALUES FOR PLOT AND CREATE PLOT 66 INPUT109 INPUT110 INPUT1L11 INPUT112 INPUT113 INPUT114 INPUT115 INPUT116 INPUT117 INPUT118 INPUT119 INPUT120 INPUT121 INPUT122 INPUT123 INPUT124 INPUT125 INPUT126 INPUT127 INPUT128 INPUT129 INPUT130 INPUT131 INPUT132 INPUT133 INPUT134 INPUT135 INPUT136 INPUT137 INPUT138 INPUT139 INPUT140 INPUT141 PLOTER PLOTER PLOTER PLOTER PLOTER PLOTER PLOTER PLOTER PLOTERIO PLOTER1L1 PLOTER12 PLOTER13 PLOTERI14 PLOTERIS PLOTERI6 PLOTERI7 PLOTERI8 PLOTERIO PLOTER20 PLOTER21 PLOTER22 PLOTER23 PLOTER24 PLOTER2S PLOTER26 PLOTER27 ODINONL WHR Cc AaAaqanNaanan e # a 200 7e0 Ao0 NF UN=2 * # # # # # SET UP MIN AND MAX LIMITS FOR PLOT FMIN(L)=FX(Lol) FMIN(2)=FX(201]) FMAX(1)=FX(191) FMAX(2)=FK(201) # ® # # # » SET UP MIN AND MAX LIMIMTS FOR HITCH AND HEAVE RATIO FMNP=F xX (2 oNWAVE) FMXP=F X (2 »NWAVE) FMNH=F x (1 oNWAVE) FMXH=F X (1 »NWAVE) 00 200 I=1,18 IF CFX(LoT) LToFMIN(L) DFMIN(L)=FX (Lol) IF (FX(L oT) .GToFMAX(1) DFMAX(L) ZF X(1 ST) IF (FX(29T) LT oFMIN(2) DFMIN(2) =Fx(2oI) IF (FX (201) GT eFMAX (2) )FMAA (2) =F X (291) IF (I1-LEeNWAVE)GO TO 200 IF (FX (oT) .LTeFMNH) FMNH3F X (151) IF (FX (191) GT eFMXH) FMXH2F X (151) IF (FX (201) .LToF MNP) FMNP2F X (291) IF (FX (291) GT eFMXP) FMXP2FX (201) CONTINUE IF (IPT.EQ.9) GO TO 800 ® # # # # COMPUTE RATIOES COL3 = (FMXH=FMNH) 7 (2.#RO) COLS% = (FMXP=FMNP) /((4,.#PI#RO) /LAMBDA) WRITE (49700) COL3,COL4 FORMAT(LH1L," HEAVE AMPLITUDE/WAVENEIGHT = "sE12.69/92Xo © % PITCH AMPLITUDE (2.#PI®wAVEHEIGHT/LAMBDA) = "sE12.6) CONTINUE NVAR(1)=10H HEAVE NVAR(2)=10H PITCH Nls2 i XO02KXA OELX = HMAx IF (NPLOT..E(),-1) CALL PLOT2(FXsFMINoF MAX »NVAR yg NFUN N15 IB,X0,D0ELX) RETURN END SUBROUTINE TRAP (F »0X»NPTS»oANS) INPUT? F ARRAY OF FUNCTIONAL VALUES OF THE INTEGRAND Ox THE X INTERVAL BETWEEN VALUES NPTS THE NUMBER OF VALUES GIVEN OUTPUT? 999 ANS THE VALUE OF THE INTEGRAL DIMENSION F(NPTS) ANS=0.0 IF (NPTS.LT.2)G0 TO 999 DO 1 I=lsNPTS ANSZANS*F (1) ANS=DX# (ANS-005#(F (1) ¢F (NPTS) )) CONTINUE RETURN END FUNCTION RMP(T»START»RISE) 67 PLOTER28 PLOTER29 PLOTER30 PLOTER31 PLOTER32 PLOTER33 PLOTER34 PLOTER3S PLOTER36 PLOTER37 PLOTER38 PLOTER39 PLOTER4O PLOTER41 PLOTER42 PLOTERS3 PLOTER44S PLOTER4S PLOTER4S6 PLOTERS7 PLOTER4S8 PLOTER49 PLOTERSO PLOTERS] PLOTERS2 PLOTERS3 PLOTERS4 PLOTERSS PLOTERS6 PLOTERS7 PLOTERSS PLOTERS9 PLOTER60 PLOTER61 PLOTER62 PLOTER63 PLOTER64 PLOTER6S PLOTER66 PLOTER67 TRAP TRAP TRAP TRAP TRAP TRAP TRAP TRAP TRAP TRAP TRAP TRAP TRAP TRAP TRAP TRAP TRAP TRAP RMP 2 3 4 5 6 7 8 9 10 ll 12 13 14 15 16 17 18 19 2 Cc Cc Cc c Cc Cc 80 99 # ® #® #@ # THIS FUNCTION IS USED TO GRADUALLY IMPLIMENT THE WAVE T CURRENT TIME START TIME TO START RAMP FRUM 0.0 TO 1.0 RISE THE LENGTH OF THE RISE FROM 0.0 TO 1.0 H=0.0 IF (T LTeSTART)GO TO 99 IF (RISE-EQ.0.0)GO TO 80 TOP=T=START Hel Ai!) IF (TOP. LT RISE) H=TOP/RISE GO JO 99 H=1. IF (T.EQoSTART)H=0.5 RMP=H RETURN END 68 RMP RMP RMP RMP RMP RMP RMP RMP RMP RMP RMP RMP RMP RMP RMP RMP RMP RMP LISTING OF COMPUTER PROGRAM FOR CALCOMP PLOTS an00 PROGRAM PLTHSP (INPUT »QUTPUT » TAPES=INPUT » TAPE6=OUTPUT 9 TAPE? 9 TAPES) ITAPE = 7 CALL CALPLT(ITAPE) STOP END SUBROUTINE CALPLT(ITAPE) DIMENSION TIME (4003) »PITCH (4003) sHEAVE (4003) * » [BUF (1000) »BWACL (4003) »CGACL (4003) . LOGICAL ACCEL 10 CAL CUMP PLOT OF PITCH AND HEAVE VERSUS TIME IREAD = S READ(IREAD.10) XAXISsYAXISP ,VAXISHohT FORMAT (8F 10.0) ACCEL = .FALSE. READ(IREAD.20) IA 9 FORMAT (110) IF (TA.EQ.1) ACCEL = .TRUE. IF (ACCEL) SCEAD(IREAD,10) YAXISbH,YAXISC CALL READT (TIME sHEAVE PITCH» BWACL sCGACL »NPTS) CALL PLOTS (IBUF 9100097) CALL PLUT(0.591.09=-3) CALL ESCALE (TIME »XAXISoNPTSo1) CALL ESCALE (HEAVE sYAXISH»NPTS91) CALL ESCALE (PITCHs YAXISP»NPTS91) IF (ACCEL) CALL ESCALE (BWACL 5 YAXISB»NPTSo1) IF (ACCEL) CALL ESCALE (CGACL 5YAXISCsNPTS91) Nl = NPTS#} N2 = NPTS¢2 N3 = NPTS¢3 CALL EAXIS(0.0,0.0,1SHTIME IN SECUNDS 9-15 9XAX1IS50.09 TIME (NL) o TIME (N2) o TIME (N3) 9HT) “CALL EAXIS(0e090ceO0o13HHEAVE IN FEET 9139 VAXISHS90009 e HEAVE (11) sHEAVE (N2) sHEAVE (N3) 9 HT) TEMP = TIME (N2) TIME (N2) = TIME (N2)/TIME(N3) HEAVE(N2) = HEAVE (N2) /REAVE (N3) CALL LINE (TIME sHEAVE »NPTS919090) TIME(N2) = TEMP XNEW = XAXTS¢3. YNEw = 1.0 CALL PLUT (xXNEW50.05=3) CALL EAAIS(00490.0,1SHTIME IN SECONDS 9=155XAXIS 90.05 e TIME (N1) o TIME (N2) 9 TIME (N3) 9HT) CALL EAXIS(00090.0e13HPITCH IN RAVo 913s YAXISP 590009 6 PITCH(NI]) sPITCH(N2) sPITCHI(N3) sAT) TIME(N2) = TIME (N2) /TIME(N3) PITCH(N2) = PITCH(N2) ZPITCH(N3) CALL LINE (TIME sPITCHsNPTS 91 9090) IF (~NOTeACCEL) GO TO 30 TIME(N2) = TEMP CALL PLOT (XNEW50.09=3) CALL EAXIS(9.0900091SHTIME IN StCUNDS9=-159XAXIS900 OsTIMECN1]) » e TIME (N2) » TIME (N3) sHT) MAIN MAIN MAIN MAIN MAIN CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALL EAXTS(2.050.0916HBUW ACCELERATION) 1659 YAXISB990.0 98WACL (N1) »CALP * BWACL (N2) sBWACL (N3) sAT) TIME (N2) = TIME (N2)/TIME(N3) -BWACL(N2) = BWACL(N2) /BWACL(N3) 69 CALP CALP CALP OBNDAMNFWNAMNFE Wh agNgNaNgN CALL LINE (TIME sBWACLoNPTS91 9000) TIME (N2) = TEMP CALL PLUT (XNEW,0.05=3) CALL EAXTS(000000091SHTIME IN SECUNDS9=159XAXIS90.009TIME (NI) » ° TIME (N2) 9 TIME (N3) HT) CALL EAXIS(0.050.091S5HCO ACCELERATION15, YAXISCs90.0,CGACL(N1) » * CGACL (N2) »CGACL (N3) sHT) 30 20 10 19 TIME (N2) = TIME (N2)/TIME(N3) CGACL(N2) = CGACL (N2) /CGACL (N3) CALL LINE (TIME »CGACL »NPTS»1 90,0) CONTINUE CALL PLOT (30.090.0,999) RETURN END SUBROUTINE READT (TIME »HEAVE »PITCHs BWACL »CGACL yNPTS) DIMENSION x (6) sHEAVE (1) sPITCH(1) o 9 TIME (1) sBWACL (1) »CGACL (1) I =0 CONTINUE I = Iel READ(9) TIME(T) 9 (X(1) 9 1=496) »BWACL (I) »CGACL (I) IF (EOF (9))109015 CONT INUE WRITE (6920) TIME(I) 9 (X (J) 9J=496) sBWACL (1) »CGACL(I) FORMAT(1H .6(F7e202X)) HEAVE(I) = xX(5) PITCH(I) = x(6) IF (1-GE-40N0) GO TO 10 GO TO S CONTINUE NPTS = l-l RETURN END SUBROUTINE EAXIS (XPAGE » YPAGE 9 IBCD »NCHAR gAXLENg ANGLE oF IRSTV 9 6 DELTAV»sDELTAUsHT) DIMENSION 1BCD(1) THIS ROUTINE wORKS LIKE THE CALCOMP AXIS WITH THE EXCEPTION THAT THE TICK MARKS ARE NOT NECCESSARILY EVERY INCH AND THE HEIGHT OF THE CHARACTERS IS INPUTTED CALL PLUT(xPAGE » YPAGE 53) ISN = ISIGN(1 »NCHAR) ISGN = SIGN(1.»DELTAV) AMIN = FIRSTV x = XPAGE Y = YPAGE XNUM = FIRSTV-DELTAV N = AXLEN/DELTAU IF (N@DELTAUCLT-AXLEN) N=Nel AMAX = AMIN*(N#DELTAV) NDIG = NDIGIT(AMINsAMAX»DELTAUsND) CONTINUE TEST = (NDIG#HT) + HT IF (TEST»GT.DELTAU) HT=HT/2. IF (TEST.GT VELTAU) GO TO 10 AYN = (1¢5#HT) BYN = (((NDIG=2) #HT) /204%.9#HT) 70 CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP CALP READ READ READ READ READ READ READ READ READ READ READ READ READ READ READ READ READ READ READ EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS tAXIS N9NANAANANAaAANA 20 30 10 20 N = N*l TANG = (90.*ANGLE) /57.22958 ANG = ANGLE/s57,2958 ST = SIN(TANG) cT = COS(TANG) S = STN(ANG) Cc = COS(ANG) 00 30 I=leN IF (I.EQ.1) GO TO 20 X = X*DELTAUSC Y = Y*DELTAU#S CALL PLOT (X9Yo2) IF (I-EQeN) GU TO 20 XT = Xe(,1#CT#ISN) YT = Yo(, 1#ST#ISN) CALL PLOT (XTsYT92) KN = X*¢AYN#CT#ISN=BYN®C YN = YeAYN#ST#ISN=BYN#®S XNUM = XNUM*¢DELTAV CALL NUMBER (XNe YNoHT »XNUM»9 ANGLE 9 ND) CALL PLOT (X9Yo3) CONTINUE XSP = (( (AXLEN/HT) 72.) —( TABS (NCHAR) /2.)) #HT YSP = 3.5#HT XT = XPAGE * XSP#C ¢ ISN#YSP#CT YT = YPAGE * XSP#S + ISN#YSP#ST CALL SYMBOL (XT YT oHT»sIBCD, ANGLE 9 LABS (NCHAR) ) RETURN END FUNCTION NDIGIT (AMIN, AMAX 9 ANUM »ND) FINDS THE NUMBER OF DIGITS NECCESSARY TO PRINT EVEN INCREMENT OF THE FUNCTION ON THE AXIS NDIGIT THE NUMBER OF PLACES IN THE ENTIRE NUMBER ND THE NUMBER OF DECIMAL PLACES ANUM THE VALUE GIVEN TO EACH INCREMENT ON THE AXIS IF (ABS (AMIN) oLT.ABS(AMAX)) GO TO 20 IF (A8S (AMIN) .EQ.ABS (AMAX) eANDeAMAX.NE.0) GO TO 20 IF (ABS (AMIN) eGT.ABS(AMAX)) GO TO 10 AMAX = 1, AMIN = -1. GO TO 20 AMAX = ABS(AMIN) IF (AMAXeLE.1.) GO TO 50 NOIV = 19 I SB i IF (AMAX/NDIV.LT.1) GO TO 40 I = I+] NDIV = NDIV#10 = wou IF (AMAX#NDIV.GT.1.) GO TO 70 I = [+l 71 EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS EAXIS NDIG NDIG NDIG NDIG NDIG NDIG NDIG NDIG NDIG NDIG NDIG NOIG NDIG NDIG NDIG NDIG NDIG NDIG NDIG NDIG NDIG NDIG NDIG NDIG NDIG NDIG NDIG NDIG NDOIG NDIG QAAAMAAANAANAIANANANN agNag0 70 a0 90 10 20 16 NDIV = NDIV#10 GO TO 60 NDIGIT = [+2 ND = I DD = FLUAT(ND) X = ANUM*(10*#DD) IX = X IF (X-FLOAT (1X) eLT220001) GO TU 90 DD = 00+1 ND = ND¢l NDIGIT = NNIGIT#l GO TO 80 CONTINUE RETURN END SUBROUTINE ESCALE (ARRAY »AXLENsNPTS» INC) FINDS THE SCALE TO BE USED ON THE AXIS = ARRAY MUST HAS THREE UNUSED POSITIONS AROAY (NPTS#1) FIRSTV ARE AY (NPTS+2) VALUES = NUMBERS) DELTAU (THE INCREMENT IN INCHES BETWEEN TICK MARKS ) AROAY (NPTS+3) DIMENSION ARRAY (1) AMIN = ARRAY (1) AMAX = ARRAY (1) ISGN = ISIGN(1>INC) INC = ITABS(INC) DO 10 Y=] eNPTSoINC IF (ACRAY (1) eLTe AMIN) AMIN=ARRAY (1) IF (ACRAY (1) eGTeAMAX) AMAX=ARRAY (I) CONTINUE AUNIT = UNIT (AMIN » AMAX o AXLENoNo ANUM) CALL AUJUST (AMIN» AMAX pAUNI T > AXLEN No ANUM) ARRAY(NPTS+1) = AMIN ARRAY (NPTS+2) = ANUM#ISGN IF (ISGN.FQe-1) ARRAY (NPTS#1) = AMAX ARRAY (NPTS+3) = AUNIT IF (ABS (ANUY) 0EQeAUNIT) ARRAY(NPTS*2) = 1.®ISGN IF (ABS(ANUM) .EQeAUNIT) ARRAY(NPTS+3) = 1. RETURN END SUBROUTINE ADJUST (AMIN» AMAX gAUNIT » AXLEN 9 Ng ANUM) GIVEN AMIN AND AMAX WHICH ARE VISTINCT VALUES» ADJUST THEM SO THAT THEY ARE EVEN MULTIPLES OF AUNIT K=1 MIN = AMIN/ANUM IF (AMIN eLT.MIN#ANUM) MIN AMIN = MIN#ANUM MAX = AMAX/ANUM IF (AMAXeGT.MAX#ANUM) MAX AMAX = MAX#ANUM TERM = AMINe (N=K) #ANUM IF (TERMeLT,AMAX) GO TO 20 MIN-1 MAX+] 72 NDIG NDIG NDIG NDIG NDIG NDIG NDIG NDIG NDIG NDIG NDIG NDIG NDIG NDIG NDIG ESCAL ESCAL ESCAL ESCAL ESCAL DELTAV (THE INCREMENT BETWEEN TICK MARKS ESCAL ESCAL ESCAL ESCAL ESCAL ESCAL ESCAL ESCAL ESCAL ESCAL ESCAL ESCAL ESCAL ESCAL ESCAL ESCAL ESCAL ESCAL ESCAL ESCAL ESCAL ESCAL ESCAL ESCAL ESCAL JUST JUST JUST JUST JUST JUST JUST JUST JUST JUST JUST JUST JUST JUST AaAaNANAaAN K = Kel GO TO 10 20 AUNIT = AXLEN/(N-Ke1) N = AXLEN/AUNI To] RETURN “END FUNCTION UNIT (AMIN, AMAX 9 AXLEN oN »ANUM) FINDS THE INCREMENT BETWEEN VALUES TO BE USED UN THE AXIS IN AS FAR AS LABELING THE TICK MARKS FINOS THE NUMBER OF DIVISIONS TO BE MADE ON THE AXIS FINDS THE SIZE IN INCHES OF THESE DIVISIUNS IF (AMINoNE.AMAX) GO TO 10 AM AM IN = AMIN-1 AX = AMAX+l 10 IF (AMAX.LT.1-AND.AMIN.GT.-1)GU TO 110 30 40 110 MIN MAX IF = AMIN = AMAX (AMAX.GToMAX) MAX=MAX¢] IF (AMINeLT.MIN) MIN=MIN-1 IF (M IF (M N IN.LT.0) NwID IN.GE.9) NwID UM = 10 MAX=MIN IF (NWID.-LT2=NUM) GO TO 60 NU GO N = M = NUM#10 To 40 NWID7 (NUM/10) MAXI ABS (MIN) IF (MINoL Te eANDeMAX.GTe0) GO TO 70 IF ( N@ (NUM/10) eL TeNwWID) N=Nel] ANUM = NUM/10. AUNIT = AXLEN/N GO TO 0 NN = IF (NN®(NUM/10) .LT.» TABS (MIN) ) NN N = IF O(N N = ANUM AUNI GO T NUM= 160 TABS (MIN) / (NUM/10) MAX (NUM/10) #(NUM/10) eLToMAX) N = Nel NeNN = NUM/10. T = AXLEN/N O 160 10 120 IF (AMAX*NUM.GT.1) GO TO 130 130 140 150 NUM GO T = NUM#10 0 120 1e/NUM AMIN®NUM AMAX#NUM IF (AMIN@®NUM.LT NI) NI=ENI=1 IF (AMAX#NUM.GT.N2) N2=N2¢1 IF (N1.NE.N2) GO TO 150 AMIN = AMIN-UNITT AMAX = AMAX=UNITT GU T9 140 N = N2-NnI ANUM = UNITT = NNel IF (AMIN. LTO. AND. AMAXoLT.0) N=N1=N2 IF (AMIN LT oOcANDeAMAX eGE VU) N=Ne=N1 AUNIT = AXLEN/N 73 JUST JUST JUST JUST JUST JUST UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT ONIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT 160 IF (NeGTeS) GO TO 170 5 N = N#2 ANUM = ANUM/2. AUNIT = AUNIT/2. GO TO 169 170 UNIT = AUNIT RETURN END 74 UNIT UNIT UNIT UNIT UNIT UNIT UNIT UNIT Copies Lh = =| = = INITIAL DISTRIBUTION WES CHONR/438 Cooper NRL 1 Code 2027 1 Code 2629 ONR/Boston ON R/Chicago ONR/Pasadena NORDA USNA 1 Tech Lib 1 Nav Sys Eng Dept 1B. Johnson 1 Bhattacheryya NAVPGSCOL 1 Library 1 T. Sarpkaya 1 J. Miller NADC NOSC 1 Library 1 Fabula 1 Hoyt NCSL/712 D. Humphreys NCEL/Code L31 NSWC, Dahlgren NUSC/Lib NAVSEA 1 SEA 0322 1 SEA 033 1 SEA 03512/Pierce 1 SEA 037 3 SEA 09G32 NAVFAC/Code 032C NAVSHIPYD PTSMH/Lib NAVSHIPYD PHILA/Lib NAVSHIPYD NORVA/Lib NAVSHIPYD CHASN/Lib NAVSHIPYD LBEACH/Lib iS) Copies 2 he — NAVSHIPYD MARE 1 Library 1 Code 250 NAVSHIPYD BREM/Lib NAVSHIPYD PEARL/Code 202.32 NAVSEC SEC 6034B SEC 6110 SEC 6114H SEC 6120 SEC 6136 SEC 6140B SEC 6144G SEC 6148 NAVSEC, NORVA/6660.03 Blount DDC AFOSR/NAM AFFOL/FYS, J. Olsen NSF/Eng Lib LC/Sci and Tech DOT/Lib TAD-491.1 MMA, Library U. of BRIDGEPORT/E. Uram U. of CAL/Dept Naval Arch, Berkeley 1 Library 1 Webster 1 Paulling 1 Wehausen U. of CAL, San Diego 1 A.T. Ellis 1 Scripps Inst Lib 1 Aero Lib 1 T.Y. Wu 1 A. Acosta CATHOLIC U. of AMER/CIVIL MECH ENG COLORADO STATE U./ENG RES CEN U. of CONNECTICUT/Scottron CORNELL U./Sears Copies FLORIDA ATLANTIC U. 1 Tech Lib 1 S. Dunne HARVARD U. 1. G. Carrier 1 Gordon McKay Lib U. of HAWAII/Bretschneider U. of ILLINOIS/J. Robertson U. of IOWA 1 Library 1 Landweber 1 Kennedy JOHN HOPKINS U./Phillips KANSAS STATE U./Nesmith U. of KANSAS/Civil Eng Lib LEHIGH U./Fritz Eng Lab Lib MIT Library Leehey Mandel Abkowitz Newman U. of MIN/ST. ANTHONY FALLS Silberman Schiebe Wetzel Song U. of MICH/NAME 1 Library 1 Ogilivie 1 Hammitt U. of NOTRE DAME 1 Eng Lib 1 Strandhagen PENN STATE/ARL/B. Parkin PRINCETON U./Mellor SIT =—= = = = Library Breslin Savitsky P.W. Brown Fridsma U. of TEXAS/ARL Lib UTAH STATE U./Jeppson aye cy ey es 76 Copies o- + = SOUTHWEST RES INST 1 Applied Mech Rev 1 Abramson STANFORD U. 1 Eng Lib 1 R. Street STANFORD RES INST/Lib U. of WASHINGTON/ARL Tech Lib WEBB INST 1 Library 1 Lewis 1 Ward WOODS HOLE/Ocean Eng WORCHESTER PIl/Tech Lib SNAME/Tech Lib BETHLEHEM STEEL/Sparrows Point BETHLEHEM STEEL/New York/Lib BOLT, BERANEK and NEWMAN/Lib GENERAL DYNAMICS, EB/Boatwright GIBBS and COX/Tech Info HYDRONAUTICS Library E. Miller A. Goodman V. Johnson 1 C.C. Hsu LOCKHEED, Sunnyvale/Waid McDONNELL DOUGLAS, Long Beach 1 J. Hess 1 T. Cebeci NEWPORT NEWS SHIPBUILDING/Lib NIELSEN ENG and RES OCEANICS ROCKWELL INTERNATIONAL/B. Ujihara SPERRY RAND/Tech Lib SUN SHIPBUILDING/Chief Naval Arch ROBERT TAGGART TRACOR —= = = = CENTER DISTRIBUTION Copies Code Name 1500 W.E. Cummins 1504 V.J. Monacella 1506 M.K. Ochi 1507 _—COD. Cieslowski 1512 J.B. Hadler 1520 RR. Wermter 1521 P. Pien 1524 = Y.T. Shen 1524 W.C. Lin 1532 G. Dobay 1532 R. Roddy 1540 W.B. Morgan 1552 = J. McCarthy 1552 N. Salvesen 1560 G. Hagen 1560 N. Hubble 1562 M. Martin 1564 J. Feldman 1568 G. Cox 1572 M.D. Ochi 1572 C.M. Lee 1576 = W.E. Smith - = =—- |} Oo - = a ee ee ey = = = ese = =| =| | = = (oe) 5214.1 Reports Distribution = 522.1 Unclassified Library (C) 1 522.2 Unclassified Library (A) Dept. Of 2 ) ne E APR 313 DTNSRDC ISSUES THREE TYPES OF REPORTS 1. DTNSRDC REPORTS, A FORMAL SERIES, CONTAIN INFORMATION OF PERMANENT TECH. NICAL VALUE. THEY CARRY A CONSECUTIVE NUMERICAL IDENTIFICATION REGARDLESS OF THEIR CLASSIFICATION OR THE ORIGINATING DEPARTMENT. 2. DEPARTMENTAL REPORTS, A SEMIFORMAL SERIES, CONTAIN INFORMATION OF A PRELIM- INARY, TEMPORARY, OR PROPRIETARY NATURE OR OF LIMITED INTEREST OR SIGNIFICANCE. THEY CARRY A DEPARTMENTAL ALPHANUMERICAL IDENTIFICATION. 3. TECHNICAL MEMORANDA, AN INFORMAL SERIES, CONTAIN TECHNICAL DOCUMENTATION OF LIMITED USE AND INTEREST. THEY ARE PRIMARILY WORKING PAPERS INTENDED FOR IN- TERNAL USE. THEY CARRY AN IDENTIFYING NUMBER WHICH INDICATES THEIR TYPE AND THE NUMERICAL CODE OF THE ORIGINATING DEPARTMENT. ANY DISTRIBUTION OUTSIDE DTNSRDC MUST BE APPROVED BY THE HEAD OF THE ORIGINATING DEPARTMENT ON A CASE-BY-CASE BASIS.