” 写作MECS 4510程序、 写作MATLAB程序设计Assignment #3ProjectMECS 4510Evolutionary ComputationLipson Pollack, Nature 406, 2000Assignment 3 3a Develop a 3D physics simulator Demonstrate a bouncing cube and a breathing cube 3b Evolve a Fixed-shape robot Manually design a robot, evolve a controller to make it move 3c Evolve both the shape and behavior of a robot Demonstrate a moving robot of various shapes and behavoors Final PresentationMock examReviewProject final presentationsHow to build a 3D simulatorBasic stepsExisting Physics SimulatorsBuild your own!Why build your own physics simulator? Simpler than learning to use other simulator Simpler to integrate with your EA Robust (wont crash) with garbage generated by an EA More versatile physics An extra line on your resumeBasic primitives Masses Point masses (no volume) Have position, velocity and acceleration Springs Massless Connect two masses Have a rest length and stiffness Apply restoring forces on the massesBasic simulator Choose a discrete time step dt Not larger than a millisecond (dt=0.001). Set global time variable T = 0 At each time step: T = T + dt Interaction step: Compute and Tally all the forces acting on each mass from springs connected it Integration step: Update position and velocit of each atom using Newtons laws of motion FсmaStick to consistent unit system, e.g. MKSRelaxation SimulationFirst: Data structures Mass Mass m (scalar) Position p (3D Vector) Velocity v (3D Vector) Acceleration a (3D Vector) External forces F (3D Vector) Spring Rest length L0 (Scalar) Spring constant k (scalar) (=EA/L for bars) Indices of two masses m1, m2Dynamic Simulation Time increment T = T +dt Interaction step For all springs: Calculate compression/extension using F=k(L-L0) For all masses: Tally forces on each mass Add gravity, collision, external forces (if any) Integration Step Calculate a = F/m (from F=ma) Update v = v + dt*a; Update p = p + v*dt Repeat forever or until reaches equilibriumdt = Simulation timestep (0.001)Dynamic = real time; geometry, mass, forces, velocities, accelerationsDynamic Simulation Time increment T = T +dt Interaction step For all springs: Calculate compression/extension using F=k(L-L0) For all masses: Tally forces on each mass Add gravity, collision, external forces (if any) Integration Step Calculate a = F/m (from F=ma) Update v = v + dt*a; If mass is fixed, set v=0 If dampening, reduce speed v=v*0.999 If colliding with ground, apply restoring force Update p = p + v*dt Repeat forever or until reaches equilibriumdt = simulation timestep (0.001)Dynamic = real time; geometry, mass, forces, velocities, accelerationsHow to handle collisions with the ground? Apply a ground reaction force FG once a node goesbelow ground (z0) FG= (0,0,-z*KG) KG is the spring coefficient of the ground, e.g. 100,000const double GRAVITY = 9.81;const double damping = 0.9999;const double DT = 0.0008; // simulation timestepconst double friction_mu_s = 1; // friction coefficient rubber-concreteconst double friction_mu_k = 0.8;// friction coefficient rubber-concreteconst double k_vertices_soft = 2000; // spring constant of the edgesconst double k_ground = 200000;// constant of the ground reaction forceDouble omega = 10;// this is arbitrary, you can use any that you see fitCreate more complex structures Add more cubes Share vertices and share springs Try making structures out of tetrahedra or other primitivesGPU- Accelerated Physis Simulations with CUDA GPU (CUDA) accelerated physics simulation sandbox capable of running complex spring/mass simulations in a fraction of the time of traditional CPU- based software while the CPU continues to perform optimizations. Capable of simulating 100,000 springs at about 4000 frames per second, or about 400,000,000 spring updates per secondint main() {Mass * m1 = sim.createMass(Vec(1, 0, 0)); // create two masses at positions (1, 0, 0), (-1, 0, 0)Mass * m2 = sim.createMass(Vec(-1, 0, 0));Spring * s1 = sim.createSpring(m1, m2); // connect them with springsCube * c1 = sim.createCube(Vec(0, 0, 3), 1); // create fully-connected cube with side length 1 and center (0, 0, 3)Plane * p = sim.createPlane(Vec(0, 0, 1), 0); // constraint plane z = 0 (constrains object above xy plane)sim.setTimeStep(0.0001) // increment simulation in 0.0001 second intervalsSim.setBreakpoint(5.0); // Ends simulation at 5 secondssim.run( ); // run simulation// perform other optimizations on the CPU while the simulation continues to run, like topology optimization, robotics, etc.return 0;}Interested in CUDA library? Contact jacob.austin@columbia.edu and Rafael Corrales Fatou rc2997@columbia.eduTitanCreate animated structures Dynamically change L0 of various springs For example, change L0=a+b*sin(ZT+c) The coefficient Z is a global constant that determines the frequencyof oscillations. Setting Z =1 will result in a period of 2S sec Chose a,b,c randomly within a reasonable range, for each spring. Keep some springs at a constant rest length by choosing b=0Simulation efficiency Measured in spring evaluations per second Per real second, not per simulated second Typical values Python: 10,000 springs/sec C++: 1.5 Million springs/sec per CPU core CUDA/GPU: 600 Million springs/sec (TitanX 2500 cores) If efficiency is low, keep structures simpleElements per second# ElementsChange material property Set spring constant to represent hard and soft materials k = 10000 Hard material k = 1000 soft material Change k as a function of spring length (nonlinear material) Change k with time? Note: Higher k Needs smaller simulation dt3D Graphics Visualization is important 3D graphics available for python, MATLAB, C++ Draw all springs and masses Draw ground plane clearly (e.g. grid in different colors) Draw external forces Use color for visualization, e.g. Change spring color with stress (e.g. red in tension, blue in compression) Change spring color depending on material type, k, etc. You dont need to draw the springs and masses ever time step For example, you can have a time-step of dt=0.0001 sec but draw the cubeonly every 100 time-stepsBoxi Xia -Solid appearance You scan make the robot look solid by shading (filling in) all theexterior triangles of the robot.General tips Choosing k and dt is trick If the structure is too wobbl k is toosmall If the structure vibrates k is too high. If the structureexplodes our timestep dt or k is too large Generall set dt to asmall step (e.g. 1E-5) and set k to a low stiffness (e.g. 10,000), andincrease them from there*. Once you get it to work, You can easily parallelize the loop in 4a and4b to use all your CPU cores. For example, in C++ use OMP. If you arereally ambitious, you can even use GPU/CUDA. Be sure to use double precision numbers for all your calculations.*Some students reported good results with k=100,000 and dt = 0.0008Validation If there is no damping or friction,energy should remain constant Plot the total energy of the cube asfunction of time (kinetic energy, potentialenergy due to gravity, potential energy inthe springs, as well as the energy relatedto the ground reaction force). The sumshould be nearly constant. If it is notconstant, you have some bugAdd friction When a mass is at or below the ground, you cansimulate friction and sliding. Assume s is acoefficient of friction and k is the kinetic coefficient Calculate the Horizontal force FH=sqrt(Fx2+Fy2) and thenormal force Fn If Fn0 (i.e. force pointing downward) then: If FH -Fn*s then zero the horizontal forces of the mass If FHt -Fn*s then reduce the horizontal force of the massby k*Fn Update the mass velocity based on the total forces,including both friction and ground reaction forcesAdd dampening Velocity dampening Multiply the velocity V by 0.999 (or similar) each time step Helps reduce vibration, wobbling, and snappingAdd dampening Velocity dampening Multiply the velocity V by 0.999 (or similar) each time step Helps reduce vibration and wobbling Add water/viscous dampening Apply drag force FD=cv2 in the opposite direction of v (c is drag coef).Anchors Anchor points by setting v=0 Constrain points to a curve or surface, e.g. set vx=0Add nonlinear behavior Replace F=k(L-L0) with some nonlinear behavior Cable F=k(L-L0) if LL0 F=0 if L= 0Final submission for 3a: Breathing cube For each spring set k and L0 = a+b*sin(Zt+c) Z is global frequency. a,b,c,k are spring parameters (k is spring coef) Cycle Period is 2S/ZAssignment 3b Parametrically evolve spring rest lengths coefficients to make amanually designed robot moveAssignment 3bEvolving Motion For each spring set k and L0 = a+b*sin(Zt+c) Z is global frequency. a,b,c,k are spring parameters (k is spring coef) Cycle period is 2S/Z Evolve locomotion Fitness: net distance travelled by center of mass Representation Direct encoding: Chromosome with a,b,c,k for each spring Indirect encoding: Chromosome specifies how a,b,c,k are determinedDirect encodinga1,b1,c1,k1a2,b2,c2,k2an,bn,cn,knFitness? Mutation? Crossover?n = number of springsIncorporating Sensors You Can change the length L0 of springs also as function of sensors Sensor types: External, e.g. light level Distance to goal Internal (proprioceptive) Actual length of spring Ground force (Whether mass is on or below ground) SpeedMore complex functions L0 = a+b*sin(Zt+c) L0 = a+b*sin(Zt+c)+d*sin(2Zt+e) L0 = a+b*sin(Zt+c)+d*sin(2Zt+e)+f*sin(3Zt+g) L0 = f(sin(Zt)) Evolve f directlyAssignment 3c Evolve morphologyDirect morphological operators Control changes Change actuation parameters a,b,c,k Morphological changes Add/remove spring Between two exiting masses Add/remove mass Change mass distribution (but keep total mass constant) Change rest length of Existing spring Some operators require cleanup and bookkeeping For example, after removing a mass you may need to remove dangling springsAnd update indices of other massesDevelopmental encodingsBase TetrahedronGrowing DesignsNervus Systems Inc.Generative Blueprints如有需要，请加QQ：99515681 或邮箱：99515681@qq.com
“
添加老师微信回复‘’官网 辅导‘’获取专业老师帮助，或点击联系老师1对1在线指导。
