Upload New File

Data_3.0/stepsize/simulation-plotting-stepsize.py

+#!/usr/bin/env python3
+# -*- coding: utf-8 -*-
+"""
+Created on Tue Apr  9 23:41:08 2019
+
+@author: administer
+"""
+import matplotlib.pyplot as plt
+import subprocess
+import os
+import glob
+import shutil
+import numpy as np
+import sys
+
+executable = "/Users/administer/Workspace/jpscore/bin/jpscore"
+
+Model = glob.glob("m*")
+for model in Model:
+        os.chdir(model)  
+        subprocess.call(["python", "makeini.py", "-f", "master_ini.xml"])     
+        os.chdir('..')
+print("inifiles are finished")
+
+# from here begins the simulation 
+
+print('from here begins the simulation ')
+#verify width list
+Model = glob.glob("m*")
+
+for model in Model:
+     os.chdir(model+'/inifiles')
+     inifiles = glob.glob("*.xml")
+     for inifile in inifiles:
+         subprocess.call([executable, "--inifile=%s"%inifile])    
+         inifiles=glob.glob('*')
+         for filename in inifiles:
+             if filename.endswith('.txt'):
+                 shutil.move(filename,'../clogginglog/')
+                 os.chdir("../clogginglog")
+                 os.rename(filename,inifile.split('size_')[1].split('.xm')[0]+'Clo'+inifile.split('d_')[1].split('_st')[0]+filename)
+                 os.chdir("../inifiles")
+  
+     os.chdir("../..")
+
+j=0
+Model = glob.glob('m*')
+
+BIG_meanNmax=[[],[],[],[],[]]   
+BIG_meanTimespan=[[],[],[],[],[]] 
+BIG_C=[[],[],[],[],[]]
+
+for model in Model:
+        
+        Nmax_w1=[]
+        Nmax_w2=[]
+        Nmax_w3=[]
+        Nmax_w4=[]
+        Nmax_w5=[]
+    
+        
+        Timespan_w1=[]
+        Timespan_w2=[]
+        Timespan_w3=[]
+        Timespan_w4=[]
+        Timespan_w5=[]
+    
+      
+        W=[0.01,0.02,0.05,0.1,0.2]
+        MeanNmax=[0,0,0,0,0,]
+        meanTimespan=[0,0,0,0,0,]
+      
+        print(model)
+        os.chdir(model+'/clogginglog')
+        Logs=glob.glob('*')
+        for log in Logs:
+            
+            
+            if log.split('Clo')[0] == '0.01':
+                if  np.loadtxt(log).shape ==(3,):
+                    Nmax_w1.append(np.loadtxt(log)[-1])
+                    MeanNmax[0]=sum(Nmax_w1)/len(Nmax_w1)
+                    Timespan_w1.append(np.loadtxt(log)[1])
+                    meanTimespan[0]=sum(Timespan_w1)/len(Timespan_w1)
+                else   :
+                    Nmax_w1.append(np.loadtxt(log)[:,2][-1])
+                    MeanNmax[0]=sum(Nmax_w1)/len(Nmax_w1)
+                    Timespan_w1.append(np.loadtxt(log)[-1][1]-np.loadtxt(log)[0][1])
+                    meanTimespan[0]=sum(Timespan_w1)/len(Timespan_w1)
+              
+                
+            if log.split('Clo')[0] == '0.02':
+                if  np.loadtxt(log).shape ==(3,):
+                    Nmax_w2.append(np.loadtxt(log)[-1])
+                    MeanNmax[1]=sum(Nmax_w2)/len(Nmax_w2)
+                    Timespan_w2.append(np.loadtxt(log)[1])
+                    meanTimespan[1]=sum(Timespan_w2)/len(Timespan_w2)
+                    
+                else  :  
+                    Nmax_w2.append(np.loadtxt(log)[:,2][-1])
+                    MeanNmax[1]=sum(Nmax_w2)/len(Nmax_w2)
+                    Timespan_w2.append(np.loadtxt(log)[-1][1]-np.loadtxt(log)[0][1])
+                    meanTimespan[1]=sum(Timespan_w2)/len(Timespan_w2)
+                    
+            if log.split('Clo')[0] == '0.05':
+                if  np.loadtxt(log).shape ==(3,):
+                    Nmax_w3.append(np.loadtxt(log)[-1])
+                    MeanNmax[2]=sum(Nmax_w3)/len(Nmax_w3)
+                    Timespan_w3.append(np.loadtxt(log)[1])
+                    meanTimespan[2]=sum(Timespan_w3)/len(Timespan_w3)
+                    
+                else   : 
+                    Nmax_w3.append(np.loadtxt(log)[:,2][-1])
+                    MeanNmax[2]=sum(Nmax_w3)/len(Nmax_w3)
+                    Timespan_w3.append(np.loadtxt(log)[-1][1]-np.loadtxt(log)[0][1])
+                    meanTimespan[2]=sum(Timespan_w3)/len(Timespan_w3)
+                    
+            if log.split('Clo')[0] == '0.1':
+                if  np.loadtxt(log).shape ==(3,):
+                    Nmax_w4.append(np.loadtxt(log)[-1])
+                    MeanNmax[3]=sum(Nmax_w4)/len(Nmax_w4)
+                    Timespan_w4.append(np.loadtxt(log)[1])
+                    meanTimespan[3]=sum(Timespan_w4)/len(Timespan_w4)
+                    
+                else :   
+                    Nmax_w4.append(np.loadtxt(log)[:,2][-1])
+                    MeanNmax[3]=sum(Nmax_w4)/len(Nmax_w4)
+                    Timespan_w4.append(np.loadtxt(log)[-1][1]-np.loadtxt(log)[0][1])
+                    meanTimespan[3]=sum(Timespan_w4)/len(Timespan_w4)
+                
+            if log.split('Clo')[0] == '0.2':      
+                if  np.loadtxt(log).shape ==(3,):
+                    Nmax_w5.append(np.loadtxt(log)[-1])
+                    MeanNmax[4]=sum(Nmax_w5)/len(Nmax_w5)
+                    Timespan_w5.append(np.loadtxt(log)[1])
+                    meanTimespan[4]=sum(Timespan_w5)/len(Timespan_w5)
+                    
+                else: 
+                    Nmax_w5.append(np.loadtxt(log)[:,2][-1])
+                    MeanNmax[4]=sum(Nmax_w5)/len(Nmax_w5)
+                    Timespan_w5.append(np.loadtxt(log)[-1][1]-np.loadtxt(log)[0][1])
+                    meanTimespan[4]=sum(Timespan_w5)/len(Timespan_w5)
+                
+              
+        print(MeanNmax)
+        print(meanTimespan)
+        BIG_meanNmax[j]=MeanNmax
+        BIG_meanTimespan[j]=meanTimespan
+        
+        C=[0,0,0,0,0]
+        for i in range(len(MeanNmax)):
+            if  meanTimespan[i]==0:
+                pass
+    
+            else:
+                C[i]=MeanNmax[i]/meanTimespan[i]
+        print(C)
+        BIG_C[j]=C
+        j=j+1
+        
+        plt.figure()
+        plt.plot(W,MeanNmax)
+        plt.savefig('W-N',dpi=300)
+        
+        plt.figure()
+        plt.grid()
+        plt.plot(W,C)
+        plt.savefig('W-Nt',dpi=300)
+        
+        
+        #plt.savefig('W-N/T')
+        
+        print('now leave'+model)
+        
+            
+        os.chdir('../..')
+        
+std1=np.std(BIG_meanNmax)
+plt.figure()
+plt.grid()
+plt.title('standard deviation is %s'%std1)
+plt.plot(W,BIG_meanNmax[0],label='M1')
+plt.plot(W,BIG_meanNmax[1],label='M0')
+plt.plot(W,BIG_meanNmax[2],label='M2')
+plt.plot(W,BIG_meanNmax[3],label='M3')
+plt.plot(W,BIG_meanNmax[4],label='M4')
+plt.xlabel('inflowrate in (m)')
+plt.ylabel('mean clogging times')
+plt.legend()
+plt.savefig('figure/multiM-stepsize-n',dpi=300)
+
+
+std2=np.std(BIG_C)
+plt.figure()
+plt.grid()
+plt.title('standard deviation is %s'%std2)
+plt.plot(W,BIG_C[0],label='M1')
+plt.plot(W,BIG_C[1],label='M0')
+plt.plot(W,BIG_C[2],label='M2')
+plt.plot(W,BIG_C[3],label='M3')
+plt.plot(W,BIG_C[4],label='M4')
+plt.xlabel('inflowrate in (m)')
+plt.ylabel('mean clogging times per seconds')
+plt.legend()
+plt.savefig('figure/multiM-stepsize-nt',dpi=300)