bin/plotForces.py

   1 #!/usr/bin/python
   2
   3 forcesfilename = 'forces/0/forces.dat'
   4
   5 import sys
   6 if len(sys.argv) != 1:
   7         print 'script assumes forces file ', forcesfilename
   8         sys.exit()
   9
  10 import re
  11 forceRegex=r"([0-9.Ee\-+]+)\s+\(+([0-9.Ee\-+]+)\s([0-9.Ee\-+]+)\s([0-9.Ee\-+]+)\)\s\(([0-9.Ee\-+]+)\s([0-9.Ee\-+]+)\s([0-9 .Ee\-+]+)\)+\s\(+([0-9.Ee\-+]+)\s([0-9.Ee\-+]+)\s([0-9.Ee\-+]+)\)\s\(([0-9.Ee\-+]+)\s([0-9.Ee\-+]+)\s([0-9.Ee\-+]+)\)+"
  12 t = []
  13 fpx = []; fpy = []; fpz = []
  14 fvx = []; fvy = []; fvz = []
  15 mpx = []; mpy = []; mpz = []
  16 mvx = []; mvy = []; mvz = []
  17
  18 pipefile=open(forcesfilename,'r')
  19 lines = pipefile.readlines()
  20
  21 for line in lines:
  22         match=re.search(forceRegex,line)
  23         if match:
  24                 t.append(float(match.group(1)))
  25                 fpx.append(float(match.group(2)))
  26                 fpy.append(float(match.group(3)))
  27                 fpz.append(float(match.group(4)))
  28                 fvx.append(float(match.group(5)))
  29                 fvy.append(float(match.group(6)))
  30                 fvz.append(float(match.group(7)))
  31                 mpx.append(float(match.group(8)))
  32                 mpy.append(float(match.group(9)))
  33                 mpz.append(float(match.group(10)))
  34                 mvx.append(float(match.group(11)))
  35                 mvy.append(float(match.group(12)))
  36                 mvz.append(float(match.group(13)))
  37
  38
  39 # combine pressure and viscous forces and moments
  40 fx = fpx
  41 fy = fpy
  42 fz = fpz
  43
  44 mx = mpx
  45 my = mpy
  46 mz = mpz
  47
  48 for i in range(1,len(t)):
  49         fx[i] += fvx[i]
  50         fy[i] += fvy[i]
  51         fz[i] += fvz[i]
  52         mx[i] += mvx[i]
  53         my[i] += mvy[i]
  54         mz[i] += mvz[i]
  55
  56 # write clean data file
  57 outForces=open('forces.dat','w')
  58
  59 for data in zip(t,fx,fy,fz):
  60         outForces.write(' '.join([str(d) for d in data])+'\n')
  61
  62 outForces.close()
  63
  64 outMoments=open('moments.dat','w')
  65
  66 for data in zip(t,mx,my,mz):
  67         outMoments.write(' '.join([str(d) for d in data])+'\n')
  68
  69 outMoments.close()
  70
  71 # plot forces
  72 import pylab
  73 pylab.xlabel('iteration')
  74 pylab.ylabel('force (N)')
  75 pylab.grid(True)
  76 #
  77 pylab.plot(t,fx,'-',label="fx")
  78 pylab.plot(t,fy,'-',label="fy")
  79 pylab.plot(t,fz,'-',label="fz")
  80
  81 #pylab.plot(t,mx,'-',label="mx")
  82 #pylab.plot(t,my,'-',label="my")
  83 #pylab.plot(t,mz,'-',label="mz")
  84
  85 pylab.legend()
  86 pylab.show()