Commit 1f69c7c0 authored by ylan's avatar ylan
Browse files

add txt files and light python chk

parent 02c20cbf
import sys
import numpy as np
#from scipy.spatial.distance import pdist # mem expansive
from scipy.spatial import Delaunay
from scipy.interpolate import interp1d
fi='ti.txt' # default input filename
fo='' # default output filename
rad0=0.91945 # radius of spheres for input
rad1=1.0 # radius of spheres for output
nkiss=2.5 # Use this kissing number to define dsep
nargv=len(sys.argv)
if nargv == 1 or nargv>4:
print("Usage: python3 %s <input filename> (no output)"%sys.argv[0])
print("Usage: python3 %s <input filename> <output_filename>"%sys.argv[0])
print("Usage: python3 %s <input filename> <output_filename> <rad0>"%sys.argv[0])
quit()
if nargv > 1:
fi=sys.argv[1]
if nargv > 2:
fo=sys.argv[2]
if nargv > 3:
rad0=sys.argv[3]
# Read input file
with open(fi,'r') as fh:
print('Input file: '+fi)
centers = np.fromfile(fh, sep=' ',dtype=float)
centers = centers.reshape((centers.size//3,3)); Nsph = centers.shape[0]
print('Nsph %d'%Nsph)
# cheap pdist
def comp_dist(P):
tri = Delaunay(P, qhull_options="Qbb C-0")
arr1 = P[tri.simplices[:,1],:] - P[tri.simplices[:,0],:]
arr2 = P[tri.simplices[:,2],:] - P[tri.simplices[:,1],:]
arr3 = P[tri.simplices[:,0],:] - P[tri.simplices[:,2],:]
arr4 = P[tri.simplices[:,3],:] - P[tri.simplices[:,0],:]
arr5 = P[tri.simplices[:,3],:] - P[tri.simplices[:,1],:]
arr6 = P[tri.simplices[:,3],:] - P[tri.simplices[:,2],:]
arr=arr1;
arr=np.append(arr,arr2,axis=0)
arr=np.append(arr,arr3,axis=0)
arr=np.append(arr,arr4,axis=0)
arr=np.append(arr,arr5,axis=0)
arr=np.append(arr,arr6,axis=0)
dist=np.sqrt(arr[:,0]**2+arr[:,1]**2+arr[:,2]**2)
return dist
dist = comp_dist(centers)
print('dtouch %f'%min(dist))
dist=np.sort(dist); Nplt=dist.shape[0];
tx = np.linspace(1,Nplt,num=Nplt,endpoint=True); tx = tx/Nsph*2
f = interp1d(tx,dist,kind='cubic')
print('dsep %f'%f(nkiss))
if 'fo' in locals() and fo != '':
# rescale domain by rad1/rad0
print('Rescale: from rad=%f to rad=%f'%(rad0,rad1))
centers = centers * rad1/rad0
# dump to output file
print('Save output into '+fo)
np.savetxt(fo,centers)
This source diff could not be displayed because it is too large. You can view the blob instead.
This source diff could not be displayed because it is too large. You can view the blob instead.
import numpy as np
from scipy.spatial.distance import pdist
rad0=0.91945
rad1=1.0
with open('pebble_centers_re2.txt','r') as fh:
centers = np.fromfile(fh, sep=' ',dtype=float)
centers = centers.reshape((centers.size//3,3))
centers = centers * rad1/rad0
np.savetxt('pebble_centers_rad_1.txt',centers)
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment