window.addEventListener('load', function(){ //JQuery function /* Thank you all for 100 likes! The simulation should work in browsers and on phones thanks to WebGL and OpenGL ES. I hope you like the simulation. Use computer for best readability of the code. Some information for those interested, these are my preferences when making 3D simulations: - Writing the JavaScript code: Sublime Text - Testing the code: Firefox browser - Making the models: Blender (and exporting as .obj and .mtl) - Finding equations of motion: MatLab */ console.log=function(){}; alert("Hi, and thanks for checking out my code!\n\nThis is a ship simulation that shows the effect of two gyroscopic stabilizers. The principle is that the moment induced by the gyroscopes cancels out the moment induced by waves.\n\nThe equations of motion were found using The Moving Frame method in dynamics with the help from Thomas Impelusso. The equations of motion were found using MatLab. Keep in mind that many simplifications have been made, and thus this is not a true simulation of what would really happen. The work presented is the result of a bachelor thesis written by me and my group.\n\nPlay around with the GUI and have fun! :-)\n\nFor those who want to learn, the maths is probably confusing, but I have tried to comment most of the code."); /* Thanks to the rest of my bachelor group: - Håkon Teigland - Andreas Flåten - Martin R. Lied - Caspar C. Smith The bachelor thesis can be read here: http://home.hib.no/prosjekter/dynamics/2017/gyroscope/ */ // ********** Setting parameters ********** const RHO = 1000; const PI = Math.PI; const g = 9.81; // Ship properties const SCALE = 27; const DRAFT = 2.578/SCALE; const LOA = 43/SCALE; const BEAM = 12/SCALE; const KEELLENGTH = 40.6/SCALE; const KEELBREADTH = 1/SCALE; const KEELHEIGHT = 0.625/SCALE; const BOTTOMBREADTH = 5.5/SCALE; const BOTTOMHEIGHT = 1.275/SCALE; // Wave properties const HEIGHT = 2.8/SCALE; const PERIOD = 5.8/Math.sqrt(SCALE); const DEPTH = 3; const LAMBDA = getWavelength(PERIOD, DEPTH); const WAVETYPE = checkDepth(HEIGHT, LAMBDA); const ZETA = HEIGHT*0.5; const WAVENUM = 2*PI/LAMBDA; const OMEGA = 2*PI*PERIOD; //Camera values const FOV = 45; const ASPECT = window.innerWidth/window.innerHeight; const NEAR = 0.1; const FAR = 2000; // ********** Creating the scene: ********** var renderer = new THREE.WebGLRenderer({ antialias: true }); //Creates a WebGL renderer to the "renderer" variable renderer.setPixelRatio( window.devicePixelRatio ); //Prevents blurry output renderer.setSize( window.innerWidth,window.innerHeight ); //Sets renderer size to the size of the window renderer.setClearColor( 0xffffff, 1); //Makes the background color of the scene white document.body.appendChild( renderer.domElement ); //Attaches renderer to DOM (initializes renderer) var scene = new THREE.Scene(); //Creates an empty scene where we are going to add our objects var camera = new THREE.PerspectiveCamera( FOV,ASPECT,NEAR,FAR ); //Creates a camera (these values have been set at the beginning) camera.position.set( 60/SCALE , 20/SCALE , 10/SCALE ); //Positions the camera in front of the ship camera.up.set( 0,0,1 ); //Rotates the camera scene.add( camera ); //Adds the camera to the scene var controls = new THREE.OrbitControls( camera, renderer.domElement ); //OrbitControls allows camera to be controlled and orbit around a target controls.minDistance = 1/SCALE; //Sets the minimum distance one can pan into the scene controls.maxDistance = 1000/SCALE; //Sets the maximum distance one can pan away from scene // ********** Initialize GUI (Graphical user interface, aka. controls) ********** var gui = new dat.GUI(); var guiParams = { //Sets the variable names and default values of the GUI RPM : 10000, gimbalAMP : 0, gimbalPeriod : 1, waves : false }; gui.add(guiParams, 'RPM', 0, 20000).name('Disk RPM'); //Adds the "Disk RPM" control to the GUI, with a range from 0-20000. gui.add(guiParams, 'gimbalAMP', 0, 90).name('Nutation AMP'); //Adds the "Nutation AMP" control to the GUI, with a range from 0-90. gui.add(guiParams, 'gimbalPeriod', 0.2, 5).name('Nutation period'); //Adds the "Nutation period" control to the GUI, with a range from 0.2-5. gui.add(guiParams, 'waves'); //Adds the "waves" button to the GUI. // ********** Loading the models ********** var ship = new THREE.Object3D(); //Creates a new threejs object named ship var mtlLoader = new THREE.MTLLoader(); //Creates a mtlLoader (to apply texture to 3d objects) mtlLoader.load( 'https://dl.dropbox.com/s/lpu2wttdnb2kuos/baat.mtl', function( materials ) //Prepare to set color to the ship { var objLoader = new THREE.OBJLoader(); //Creates a object loader (to load the ship) objLoader.setMaterials( materials ); //Sets color to the ship objLoader.load( 'https://dl.dropbox.com/s/o6s45mo3r1d8juf/baat.obj', function ( object )//Loads ship { ship.add( object ); //Adds ship to scene }); }); //Same approach for all objects var gyroCribPS = new THREE.Object3D(); mtlLoader = new THREE.MTLLoader(); mtlLoader.load( 'https://dl.dropbox.com/s/m4gunv3p1iktf1n/crib.mtl', function( materials ) { var objLoader = new THREE.OBJLoader(); objLoader.setMaterials( materials ); objLoader.load( 'https://dl.dropbox.com/s/27zgge6modbv7qn/crib.obj', function ( object ) { gyroCribPS.add( object ); }); }); var gyroCribSB = new THREE.Object3D(); mtlLoader = new THREE.MTLLoader(); mtlLoader.load( 'https://dl.dropbox.com/s/m4gunv3p1iktf1n/crib.mtl', function( materials ) { var objLoader = new THREE.OBJLoader(); objLoader.setMaterials( materials ); objLoader.load( 'https://dl.dropbox.com/s/27zgge6modbv7qn/crib.obj', function ( object ) { gyroCribSB.add( object ); }); }); var pinionPS = new THREE.Object3D(); mtlLoader = new THREE.MTLLoader(); mtlLoader.load( 'https://dl.dropbox.com/s/xrhrp88g2ohsfj4/pinion.mtl', function( materials ) { var objLoader = new THREE.OBJLoader(); objLoader.setMaterials( materials ); objLoader.load( 'https://dl.dropbox.com/s/fa1413qjqeb4ybm/pinion.obj', function ( object ) { pinionPS.add( object ); }); }); var pinionSB = new THREE.Object3D(); mtlLoader = new THREE.MTLLoader(); mtlLoader.load( 'https://dl.dropbox.com/s/xrhrp88g2ohsfj4/pinion.mtl', function( materials ) { var objLoader = new THREE.OBJLoader(); objLoader.setMaterials( materials ); objLoader.load( 'https://dl.dropbox.com/s/fa1413qjqeb4ybm/pinion.obj', function ( object ) { pinionSB.add( object ); }); }); var gimbalPS = new THREE.Object3D(); mtlLoader = new THREE.MTLLoader(); mtlLoader.load( 'https://dl.dropbox.com/s/2eycbpdf1y0c6ih/gimbal.mtl', function( materials ) { var objLoader = new THREE.OBJLoader(); objLoader.setMaterials( materials ); objLoader.load( 'https://dl.dropbox.com/s/espi9yqxnpybt8o/gimbal.obj', function ( object ) { gimbalPS.add( object ); }); }); var gimbalSB = new THREE.Object3D(); mtlLoader = new THREE.MTLLoader(); mtlLoader.load( 'https://dl.dropbox.com/s/2eycbpdf1y0c6ih/gimbal.mtl', function( materials ) { var objLoader = new THREE.OBJLoader(); objLoader.setMaterials( materials ); objLoader.load( 'https://dl.dropbox.com/s/espi9yqxnpybt8o/gimbal.obj', function ( object ) { gimbalSB.add( object ); }); }); var flywheelPS = new THREE.Object3D(); mtlLoader.load( 'https://dl.dropbox.com/s/962oi7hebnmfp54/flywheel.mtl', function( materials ) { var objLoader = new THREE.OBJLoader(); objLoader.setMaterials( materials ); objLoader.load( 'https://dl.dropbox.com/s/ie6bzfmzr02yx5r/flywheel.obj', function ( object ) { flywheelPS.add( object ); }); }); var flywheelSB = new THREE.Object3D(); mtlLoader.load( 'https://dl.dropbox.com/s/962oi7hebnmfp54/flywheel.mtl', function( materials ) { var objLoader = new THREE.OBJLoader(); objLoader.setMaterials( materials ); objLoader.load( 'https://dl.dropbox.com/s/ie6bzfmzr02yx5r/flywheel.obj', function ( object ) { flywheelSB.add( object ); }); }); // ********** Setting the positions ********** gyroCribPS.position.y = 2/SCALE; //Sets the y-position of the crib gyroCribPS.position.z = 5.223/SCALE; //Sets the z-position of the crib gyroCribPS.rotation.z = PI; //Rotates one of the gimbals 180degrees so the gyros nutate opposite directions gyroCribSB.position.y = -gyroCribPS.position.y; //Cribs are symmetric about the axes gyroCribSB.position.z = gyroCribPS.position.z; gimbalPS.position.z = 2.16/SCALE; //Sets the z-position of the gimbals gimbalSB.position.z = gimbalPS.position.z; pinionPS.position.z = 0.795/SCALE; //Sets the z-position of the pinions pinionSB.position.z = pinionPS.position.z; var COF = new THREE.Object3D(); //Adding a center of focus variable to spawn the ship ship.position.z = -DRAFT; COF.add( ship ); scene.add( COF ); //Adds the ship to the scene ship.add( gyroCribPS ); //Adds the crib to the ship ship.add( gyroCribSB ); gyroCribPS.add( gimbalPS ); //Adds the gimbal to the crib gyroCribSB.add( gimbalSB ); gyroCribPS.add( pinionPS ); //Adds the pinion to the crib gyroCribSB.add( pinionSB ); gimbalPS.add( flywheelPS ); //Adds the flywheel to the gimbal gimbalSB.add( flywheelSB ); // We don't have to adjust the y-value of the gimbals and pinions, // because they are part of the crib, and their y-value is 0 relative to the crib. // ********** Making the x,y,z- axis ********** var material = new THREE.LineBasicMaterial({ //ThreeJS function for setting line color color: 0x0000ff //Setting line-color to blue }); var linegeo = new THREE.Geometry(); //Setting the geometry of the lines (X,Y,Z) linegeo.vertices.push( new THREE.Vector3( 25/SCALE, 0, 0 ), new THREE.Vector3( 0, 0, 0 ), new THREE.Vector3( 0, 25/SCALE, 0 ), //Lines are drawn between each vertices new THREE.Vector3( 0, 0, 0 ), new THREE.Vector3( 0, 0, 25/SCALE ) ); var line = new THREE.Line( linegeo, material ); //Drawing the line and setting the color scene.add( line ); //Adding the line to our scene. // ********** Adding light ********** var ambientLight = new THREE.AmbientLight(0xffffff, 0.4); scene.add(ambientLight); //Adding ambient light var light1 = new THREE.DirectionalLight( 0xffffbb, 1 ); light1.position.set( 1, 1, 1 ); scene.add( light1 ); //Adding directional light // ********** Waves ********** var geometry, points = []; var width = 100/SCALE; //width of water var height = 200/SCALE; //length of water var widthSegments = 200; //Needed for waveish plane var heightSegments = 200; geometry = new THREE.PlaneGeometry( width, height, widthSegments, heightSegments ); //ThreeJS function to create plane geometry var oceanmat = new THREE.MeshLambertMaterial( { //Sets color and material attributes for waves color: 0x2ECCFA, opacity: 0.5, transparent: true, side: THREE.DoubleSide } ); var ocean = new THREE.Mesh( geometry, oceanmat ); //Creates a mesh containing of the geometry and oceanmaterial just added scene.add( ocean ); //Adds ocean to scene // ********** Dynamics ********** var t = 0; //Time var dt = 0.01; //Timestep var w = [0,0,0]; var q = [0,0,0]; var phi = 0; var psi = 0; var rotMat = eyeMat(3,3); //Calls eyeMat function, creates an identity matrix var shipData, roll, pitch, yaw, rolld, rollvel, dpsi, RPM, phiPin; var render = function () { requestAnimationFrame(render); //render function is called every frame RPM = guiParams.RPM; //Sets the value of RPM to the GUI value dpsi = (RPM/60)*2*PI; //Calculates angular velocity of disk psi += dpsi*dt; //New angular placement of disk = Old placement + psi shipData = shipMotions(q,rotMat,t,dt); //Calls the shipMotions function and stores return data in the shipData variable //shipMotions returns rotMat: reconstructed rotation matrix, w: angular velocities, shipAngle: ship angles and phi: nutation angle w = shipData.w; rotMat = shipData.rotMat; roll = shipData.shipAngle[0]; pitch = shipData.shipAngle[1]; yaw = shipData.shipAngle[2]; phi = shipData.phi; phiPin = -0.6*1.85*phi; //Calculates how much pinions should rotate rolld = roll*180/PI; //Converts roll angle from rad to degree rollvel = w[0]*180/PI; //Converts roll velocity from rad/s to degree/s COF.rotation.x = roll; //Rotate the COF object (which contains the ship) COF.rotation.y = pitch; COF.rotation.z = yaw; gimbalPS.rotation.y = phi; //Rotates gimbals by amount: phi gimbalSB.rotation.y = phi; pinionPS.rotation.y = phiPin; //Rotates pinions by amount: phiPin pinionSB.rotation.y = phiPin; flywheelPS.rotation.z = psi; //Rotates flywheels by amount: psi flywheelSB.rotation.z = psi; q = [w[0], w[1], w[2]]; //Stores angular velocity values in "q". t += dt; //Increases time updatePoints(t); //Calls the 'updatePoints' function to update waves with t as input renderer.render(scene, camera); //We need this in the loop to actually perform the rendering }; render(); //Calls the render function again and again function updatePoints(t) { //Update waves for(var i = 0; i < ocean.geometry.vertices.length; i++) { var v = ocean.geometry.vertices[i]; var y = v.y; if ( guiParams.waves ) { // If waves are enabled: v.z = waveProfile(ZETA, OMEGA, WAVENUM, t, y); //Calls the waveProfile function } else { v.z = 0; } ocean.geometry.verticesNeedUpdate = true; } } function shipMotions(q,rotMat,t,dt) { //ShipMotions-function with four inputs q = RK4(t,q,dt); //Calls the RK4 (Runga Kutta) function. RK4 returns a 1x3 array containing w11,w12,w13 and stores it in q. var w = [ q[0], q[1], q[2] ]; //Stores the angular velocities, w11, w12 and w13 rotMat = getRotMat(rotMat,w,dt); //Reconstruction of rotation matrix var shipAngle = rotMatToAngle(rotMat); //Transforms reconstructed rotation matrix to angle var gimbalData = getNutation(t); //Calls getNutation function to get phi var phi = gimbalData.phi; //Stores phi value in phi return { //Returns reconstructed rotation matrix, angular velocities, shipAngle and phi. rotMat : rotMat, w : w, shipAngle : shipAngle, phi : phi }; } function RK4(x,y,h) //Runge Kutta function, where inputs: x=t, y=q, h=dt { var arr = [], k = [], l = [], m = [], n = []; for (var i = 0; i < y.length; i++) { k[i] = h*myFunction( i, x, y); //Calculates h*k1, h*m1, h*n1 } for ( i = 0; i < y.length; i++) { arr = constMatMult(0.5,k); // Multiplies h*k1, h*m1 and h*n1 by 0.5 l[i] = h*myFunction( i, x+(h/2), MatMatAdd(arr,y)); //Calculates h*k2, h*m2, h*n2. MatMatAdd = (h*k1*0.5)+y } for ( i = 0; i < y.length; i++) { arr = constMatMult(0.5,l); // Multiplies h*k2, h*m2 and h*n2 by 0.5 m[i] = h*myFunction( i, x+(h/2), MatMatAdd(arr,y)); //Calculates h*k3, h*m3, h*n3. MatMatAdd = (h*k2*0.5)+y } for ( i = 0; i < y.length; i++) { n[i] = h*myFunction( i, x+h, MatMatAdd(m,y)); //Calculates h*k4, h*m4, h*n4. MatMatAdd = (h*k3)+y } for ( i = 0; i < y.length; i++) { y[i] = y[i] + (k[i] + 2*l[i] + 2*m[i] + n[i])/6; //Calculates w11,w12,w13 } return y; //returns the angular velocity values [w11,w12,w13] } function myFunction(flag, t, dq) { var w11 = dq[0]; //Ship roll velocity var w12 = dq[1]; //Ship pitch velocity var w13 = dq[2]; //Ship yaw velocity const m1 = 318000/SCALE/SCALE/SCALE; const m4 = 0.0341; const m5 = 0.1453; const DELTA = m1 + 2*m4 + 2*m5; const J11 = 0.649579; const J12 = 3.549014; //Ship inertias const J13 = 3.590727; const J42 = 100.32*Math.pow(10,-6); //Gimbal inertia const J53 = 74.44*Math.pow(10,-6); // Disk inertia //We must calculate phi4, dphi4 and dpsi5 inside this function, because //they are used to calculate the equations of motion. var phi4Data = getNutation(t); //Calls the getNutation function for the gimbal with return values "phi" and "dphi" var phi4 = phi4Data.phi; //Stores the phi value obtained from getNutation var dphi4 = phi4Data.dphi; //Stores the dphi value obtained from getNutation var RPM = guiParams.RPM; //Sets the RPM value to the GUI value var dpsi5 = RPM/60*2*PI; //Calculates angular velocity of the disk // Wave forces and moments var Fw1 = 0; var Fw2 = 0; var Fw3 = 0; var Mw1 = 0; var Mw2 = 0; var Mw3 = 0; if ( guiParams.waves ) { //If waves are enabled, Mw1 get this value: Mw1 = 100 * ZETA * Math.sin(2 * PI / PERIOD * t); } // Forces var Fg = -DELTA*g; var Fb = -Fg; var F2 = Fw2; var F3 = Fg + Fb + Fw3; // Rightning Moment var GMT = 9.85/SCALE; var GML = 80.513/SCALE; var Mb1 = Fb*calculateRightingArm(GMT, COF.rotation.x); var Mb2 = -Fb*calculateRightingArm(GML, COF.rotation.y); // Damping var Broll = 2; var Bpitch = 3; var Byaw = 3; var Md1 = -Broll*w11; var Md2 = -Bpitch*w12; var Md3 = -Byaw*w13; // Total moments var M1 = Mw1 + Mb1 + Md1; var M2 = Mw2 + Mb2 + Md2; var M3 = Mw3 + Md3; //Equations obtained in MatLab! var ddq; if (flag == 0) { // dw1dt equation (to calculate k1,k2,k3,k4, w11) ddq = M1; ddq -= J13*w12*w13 - J12*w12*w13 + J53*w12*w13*Math.cos(phi4)*Math.cos(phi4) + 2*J53*Math.cos(phi4)*dphi4*dpsi5 + J53*w11*Math.sin(2*phi4)*dphi4; ddq /= J11 + 2*J42 + J53 + J53*Math.sin(phi4)*Math.sin(phi4); } if (flag == 1) { // dw2dt equation (to calculate n1,n2,n3,n4, w12) ddq = M2; ddq -= J11*w11*w13 - J13*w11*w13 - J53*w11*w13 + 2*J53*w11*w13*Math.sin(phi4)*Math.sin(phi4) + 2*J53*w13*Math.sin(phi4)*dpsi5; ddq /= J12 + 2*J42 + J53; } if (flag == 2) { // dw3dt equation (to calculate m1,m2,m3,m4, w13) ddq = M3; ddq -= J12*w11*w12 - J11*w11*w12 - J53*w11*w12 + J53*w11*w12*Math.cos(phi4)*Math.cos(phi4) - 2*J53*w12*Math.sin(phi4)*dpsi5 - 2*J53*w13*Math.cos(phi4)*Math.sin(phi4)*dphi4; ddq /= J13 + 2*J42 + 2*J53 - J53*Math.sin(phi4)*Math.sin(phi4); } return ddq; // Returns this value when RK4 calls "MyFunction" } function getNutation( t ) { var Ag = (PI/180)*guiParams.gimbalAMP; //Converts the gimbal amplitude to rads/s var Tg = guiParams.gimbalPeriod; //Stores the gimbal period var phi = Ag*Math.sin(2*(Math.PI/Tg)*t); var dphi = 2*(PI/Tg)*Ag*Math.cos(2*(PI/Tg)*t); //Calculating phi and dphi (equations from Oscar's paper) return { phi: phi, dphi: dphi //Returning phi and dphi values }; } //Below are a few matrix managing functions function zerosMat(dimensions) { var mat = []; for (var i = 0; i < dimensions[0]; ++i) { mat.push(dimensions.length == 1 ? 0 : zerosMat(dimensions.slice(1))); } return mat; } function eyeMat(rows,cols) { var Mat = zerosMat([rows,cols]); for ( var i=0; i 0.5) { return deep; } var len = deep; var temp = len; var diff = len; while ( Math.abs(diff) > 0.0001) { temp = T*T*tanh(2*PI*h/len)*g/(2*PI); diff = len - temp; len = temp; } return len; } function checkDepth( h, lambda ) { if ( h/lambda > 0.5) { return 'd'; } if ( h/lambda < 0.05 ) { return 's'; } return 'f'; } function waveProfile( zeta, omega, k, t, x ) { return zeta*Math.sin(omega*t - k*x); } function dynamicPressure( type, h, zeta, omega, k, t, x, z) { var Pd = RHO*g*zeta*Math.sin(omega*t - k*x); if (type == 'd') { Pd = Pd*Math.exp(k*z); } else if (type == 'f') { Pd = Pd*cosh(k*(z + h))/cosh(k*h); } return Pd; } function calculateRightingArm(GM, angle) { return GM*Math.sin(angle); } //Some math function tanh(Angle) { // Angle in radians var s=sinh(Angle); var c=cosh(Angle); return s/c; } function sinh(Angle) { // Angle in radians var e=Math.E; var p=Math.pow(e,Angle); var n=1/p; return (p-n)/2; } function cosh(Angle) { // Angle in radians var e=Math.E; var p=Math.pow(e,Angle); var n=1/p; return (p*1+n*1)/2; } }); //Anyone got this far?;-)

⛴3D Gyroscopic stabilizers⛴