tutorials/mesh/moveDynamicMesh/simpleHarmonicMotion/shm.gnuplot

   1 #------------------------------------------------------------------------------
   2 # =========                 |
   3 # \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4 #  \\    /   O peration     |
   5 #   \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
   6 #    \\/     M anipulation  |
   7 #-------------------------------------------------------------------------------
   8 # License
   9 #     This file is part of OpenFOAM.
  10 #
  11 #     OpenFOAM is free software: you can redistribute it and/or modify it
  12 #     under the terms of the GNU General Public License as published by
  13 #     the Free Software Foundation, either version 3 of the License, or
  14 #     (at your option) any later version.
  15 #
  16 #     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
  17 #     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
  18 #     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
  19 #     for more details.
  20 #
  21 #     You should have received a copy of the GNU General Public License
  22 #     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
  23 #
  24 # Script
  25 #     shm.gnuplot
  26 #
  27 # Description
  28 #     Creates an .eps graph of OpenFOAM results vs analytical solution
  29 #     for a simple harmonic motion dynamicMesh case
  30 #
  31 #------------------------------------------------------------------------------
  32
  33 reset
  34
  35 set samples 2000
  36
  37 k = 5000.0
  38 m = 9.6
  39 c = 50.0
  40 a = -0.1
  41
  42 omega = sqrt(k/m)
  43 zeta = c/(2.0*m*omega)
  44
  45 phi = atan((sqrt(1.0 - zeta**2))/zeta)
  46 A = a/sin(phi)
  47
  48 pos(A, t, omega, phi, zeta) = A*exp(-zeta*omega*t)*sin(sqrt(1-zeta**2)*omega*t + phi)
  49 vel(A, t, omega, phi, zeta) = \
  50 A*exp(-zeta*omega*t)*\
  51 ( \
  52   sqrt(1-zeta**2)*omega*cos(sqrt(1-zeta**2)*omega*t + phi) \
  53 - zeta*omega*sin(sqrt(1-zeta**2)*omega*t + phi) \
  54 )
  55
  56 set xlabel "Time/[s]"
  57 set ylabel "Position"
  58
  59 set ytics nomirror
  60 set y2tics
  61
  62 set yrange [-0.1:0.1]
  63 set y2range [-2:2]
  64
  65 set xzeroaxis
  66
  67 set terminal postscript eps color enhanced solid
  68 set output "shm.eps"
  69
  70 plot \
  71     "t_vs_cm" u 1:($2 - 0.6) w l t "Simulation, centre of mass relative to start", \
  72     pos(A, x, omega, phi, zeta) w l t "Analytical solution, centre of mass", \
  73     "t_vs_lv" u 1:2 w l axes x1y2 t "Simulation, vertical velocity", \
  74     vel(A, x, omega, phi, zeta) w l axes x1y2 t "Analytical solution, vertical velocity"
  75
  76 #------------------------------------------------------------------------------