From 620a42e0e8ebb2b617edab1e5433073f610fa047 Mon Sep 17 00:00:00 2001 From: David Date: Sat, 19 Sep 2020 21:39:09 +0200 Subject: [PATCH 1/2] Merge changes from BOUT-dev --- boututils/analyse_equil_2.py | 12 ++-- boututils/bunch.py | 6 ++ boututils/check_scaling.py | 11 ++-- boututils/datafile.py | 7 ++- boututils/efit_analyzer.py | 2 +- boututils/geqdsk.py | 2 +- boututils/moment_xyzt.py | 16 ++--- boututils/read_geqdsk.py | 6 +- boututils/run_wrapper.py | 117 +++++++++++++++++++++++++---------- boututils/surface_average.py | 32 ++++------ boututils/volume_integral.py | 21 +++---- 11 files changed, 133 insertions(+), 99 deletions(-) create mode 100644 boututils/bunch.py diff --git a/boututils/analyse_equil_2.py b/boututils/analyse_equil_2.py index 93b98fc..d315c6a 100644 --- a/boututils/analyse_equil_2.py +++ b/boututils/analyse_equil_2.py @@ -12,7 +12,6 @@ from past.utils import old_div import numpy -from bunch import Bunch from . import local_min_max from scipy.interpolate import RectBivariateSpline from matplotlib.pyplot import contour, gradient, annotate, plot, draw @@ -35,8 +34,8 @@ def analyse_equil(F, R, Z): Returns ------- - bunch - A structure of critical points containing: + object + An object of critical points containing: n_opoint, n_xpoint - Number of O- and X-points primary_opt - Index of plasma centre O-point @@ -169,8 +168,7 @@ def analyse_equil(F, R, Z): print("Number of O-points: "+str(n_opoint)) if n_opoint == 0 : - print("No O-points! Giving up on this equilibrium") - return Bunch(n_opoint=0, n_xpoint=0, primary_opt=-1) + raise RuntimeError("No O-points! Giving up on this equilibrium") #;;;;;;;;;;;;;; Find plasma centre ;;;;;;;;;;;;;;;;;;; @@ -259,8 +257,8 @@ def analyse_equil(F, R, Z): - #;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; - # Put results into a structure + #;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; + # Put results into a structure result = Bunch(n_opoint=n_opoint, n_xpoint=n_xpoint, # Number of O- and X-points primary_opt=primary_opt, # Which O-point is the plasma centre diff --git a/boututils/bunch.py b/boututils/bunch.py new file mode 100644 index 0000000..2bc1ca0 --- /dev/null +++ b/boututils/bunch.py @@ -0,0 +1,6 @@ +# what we need from bunch + +class Bunch: + def __init__(self, **dict): + for k in dict: + setattr(self, k, dict[k]) diff --git a/boututils/check_scaling.py b/boututils/check_scaling.py index be22fc0..af59b0b 100644 --- a/boututils/check_scaling.py +++ b/boututils/check_scaling.py @@ -44,15 +44,14 @@ def check_order(error_list, expected_order, tolerance=2.e-1, spacing=None): if len(error_list) < 2: raise RuntimeError("Expected at least 2 data points to calculate error") - success=True + success = True + for i in range(len(error_list)-1): - if spacing is None: - actual_order = log(errors[i] / errors[i+1]) / log(2) - else: - actual_order = log(errors[i] / errors[i+1]) / log(grid_spacing[i] / grid_spacing[i+1]) + grid_spacing = 2 if spacing is None else spacing[i] / spacing[i+1] + actual_order = log(error_list[i] / error_list[i+1]) / log(grid_spacing) if not isclose(actual_order, expected_order, atol=tolerance, rtol=0): - success=False + success = False return success diff --git a/boututils/datafile.py b/boututils/datafile.py index e8d385b..00bb987 100644 --- a/boututils/datafile.py +++ b/boututils/datafile.py @@ -729,11 +729,12 @@ def read(self, name, ranges=None, asBoutArray=True): if var is None: return None - if asBoutArray: - attributes = self.attributes(n) + attributes = self.attributes(n) if asBoutArray else {} + + time_dependent = attributes.get("bout_type", "none").endswith("_t") ndims = len(var.shape) - if ndims == 1 and var.shape[0] == 1: + if ndims == 1 and var.shape[0] == 1 and not time_dependent: data = var if asBoutArray: data = BoutArray(data, attributes=attributes) diff --git a/boututils/efit_analyzer.py b/boututils/efit_analyzer.py index 0b15f42..6577b78 100644 --- a/boututils/efit_analyzer.py +++ b/boututils/efit_analyzer.py @@ -4,7 +4,7 @@ from builtins import range from past.utils import old_div import numpy as np -from bunch import Bunch +from boututils.bunch import Bunch from .radial_grid import radial_grid from .analyse_equil_2 import analyse_equil from pylab import (cm, clabel, contour, draw, legend, plot, setp, show, diff --git a/boututils/geqdsk.py b/boututils/geqdsk.py index 4335b03..d63caf1 100755 --- a/boututils/geqdsk.py +++ b/boututils/geqdsk.py @@ -1,4 +1,4 @@ -#!/usr/bin/env python +#!/usr/bin/env python3 from __future__ import print_function from __future__ import division diff --git a/boututils/moment_xyzt.py b/boututils/moment_xyzt.py index 8d7dc4c..675d1f1 100644 --- a/boututils/moment_xyzt.py +++ b/boututils/moment_xyzt.py @@ -3,8 +3,7 @@ from builtins import range from past.utils import old_div import numpy as np -from bunch import Bunch - +from boututils.bunch import Bunch def RMSvalue( vec1d): #; @@ -16,9 +15,9 @@ def RMSvalue( vec1d): valrms=np.sqrt(old_div(np.sum((vec1d-valav)**2),nel)) acvec=vec1d-valav - - return Bunch(valrms=valrms, valav=valav, acvec=acvec) - + return Bunch(valrms=valrms, + valav=valav, + acvec=acvec) @@ -30,8 +29,6 @@ def moment_xyzt( sig_xyzt, *args):#rms=None, dc=None, ac=None): #; -AC (DC subtracted out), i.e., a function of (x,y,z,t) #;------------------------------------------------------------------- - try: # return to caller - d = np.shape(sig_xyzt) if np.size(d) != 4 : print("Error: Variable must be 4D (x,y,z,t)") @@ -73,8 +70,5 @@ def moment_xyzt( sig_xyzt, *args):#rms=None, dc=None, ac=None): res.ac = ac if 'RMS' not in args and 'DC' not in args and 'AC' not in args : - print('Wrong argument') + raise RuntimeError('Wrong argument') return res - except: - print('moment_xyz failed') - return diff --git a/boututils/read_geqdsk.py b/boututils/read_geqdsk.py index 755b3f1..f665f59 100644 --- a/boututils/read_geqdsk.py +++ b/boututils/read_geqdsk.py @@ -1,8 +1,8 @@ from __future__ import print_function from builtins import range import numpy -from bunch import Bunch from geqdsk import Geqdsk +from boututils.bunch import Bunch def read_geqdsk (file): @@ -67,8 +67,6 @@ def read_geqdsk (file): for j in range(0,nyefit): r[i,j] = rgrid1 + xdim*i/(nxefit-1) z[i,j] = (zmid-0.5*zdim) + zdim*j/(nyefit-1) - - f=f.T @@ -89,4 +87,4 @@ def read_geqdsk (file): pres=pres, # Plasma pressure in nt/m^2 on uniform flux grid qpsi=qpsi, # q values on uniform flux grid nbdry=nbdry, rbdry=rbdry, zbdry=zbdry, # Plasma boundary - nlim=nlim, xlim=xlim, ylim=ylim) # Wall boundary \ No newline at end of file + nlim=nlim, xlim=xlim, ylim=ylim) # Wall boundary diff --git a/boututils/run_wrapper.py b/boututils/run_wrapper.py index c25e343..9a91d2c 100644 --- a/boututils/run_wrapper.py +++ b/boututils/run_wrapper.py @@ -2,21 +2,20 @@ from builtins import str import os +import pathlib import re import subprocess +from subprocess import call, Popen, STDOUT, PIPE -try: - # Python 2.4 onwards - from subprocess import call, Popen, STDOUT, PIPE - lib = "call" -except ImportError: - # FIXME: drop support for python < 2.4! - # Use os.system (depreciated) - from os import popen4, system - lib = "system" +if os.name == "nt": + # Default on Windows + DEFAULT_MPIRUN = "mpiexec.exe -n" +else: + DEFAULT_MPIRUN = "mpirun -np" -def getmpirun(default="mpirun -np"): + +def getmpirun(default=DEFAULT_MPIRUN): """Return environment variable named MPIRUN, if it exists else return a default mpirun command @@ -53,27 +52,20 @@ def shell(command, pipe=False): """ output = None status = 0 - if lib == "system": - if pipe: - handle = popen4(command) - output = handle[1].read() - else: - status = system(command) + if pipe: + child = Popen(command, stderr=STDOUT, stdout=PIPE, shell=True) + # This returns a b'string' which is casted to string in + # python 2. However, as we want to use f.write() in our + # runtest, we cast this to utf-8 here + output = child.stdout.read().decode("utf-8", "ignore") + # Wait for the process to finish. Note that child.wait() + # would have deadlocked the system as stdout is PIPEd, we + # therefore use communicate, which in the end also waits for + # the process to finish + child.communicate() + status = child.returncode else: - if pipe: - child = Popen(command, stderr=STDOUT, stdout=PIPE, shell=True) - # This returns a b'string' which is casted to string in - # python 2. However, as we want to use f.write() in our - # runtest, we cast this to utf-8 here - output = child.stdout.read().decode("utf-8", "ignore") - # Wait for the process to finish. Note that child.wait() - # would have deadlocked the system as stdout is PIPEd, we - # therefore use communicate, which in the end also waits for - # the process to finish - child.communicate() - status = child.returncode - else: - status = call(command, shell=True) + status = call(command, shell=True) return status, output @@ -93,6 +85,18 @@ def determineNumberOfCPUs(): The number of CPUs """ + # cpuset + # cpuset may restrict the number of *available* processors + try: + m = re.search(r'(?m)^Cpus_allowed:\s*(.*)$', + open('/proc/self/status').read()) + if m: + res = bin(int(m.group(1).replace(',', ''), 16)).count('1') + if res > 0: + return res + except IOError: + pass + # Python 2.6+ try: import multiprocessing @@ -228,8 +232,12 @@ def launch(command, runcmd=None, nproc=None, mthread=None, cmd = cmd + " > "+output if mthread is not None: - cmd = "OMP_NUM_THREADS={j} ".format(j=mthread)+cmd - + if os.name == "nt": + # We're on windows, so we have to do it a little different + cmd = 'cmd /C "set OMP_NUM_THREADS={} && {}"'.format(mthread, cmd) + else: + cmd = "OMP_NUM_THREADS={} {}".format(mthread, cmd) + if verbose == True: print(cmd) @@ -277,3 +285,48 @@ def launch_safe(command, *args, **kwargs): "Output was\n\n%s"% (s,command,out)) return s, out + + +def build_and_log(test): + """Run make and redirect the output to a log file. Prints input + + On Windows, does nothing because executable should have already + been built + + """ + + if os.name == "nt": + return + + print("Making {}".format(test)) + + if os.path.exists("makefile") or os.path.exists("Makefile"): + return shell_safe("make > make.log") + + ctest_filename = "CTestTestfile.cmake" + if not os.path.exists(ctest_filename): + raise RuntimeError("Could not build: no makefile and no CMake files detected") + + # We're using CMake, but we need to know the target name. If + # bout_add_integrated_test was used (which it should have been!), + # then the test name is the same as the target name + with open(ctest_filename, "r") as f: + contents = f.read() + match = re.search("add_test.(.*) ", contents) + if match is None: + raise RuntimeError("Using CMake, but could not determine test name") + test_name = match.group(1) + + # Now we need to find the build directory. It'll be the first + # parent containing CMakeCache.txt + here = pathlib.Path(".").absolute() + for parent in here.parents: + if (parent / "CMakeCache.txt").exists(): + return shell_safe( + "cmake --build {} --target {} > make.log".format(parent, test_name) + ) + + # We've just looked up the entire directory structure and not + # found the build directory, this could happen if CMakeCache was + # deleted, in which case we can't build anyway + raise RuntimeError("Using CMake, but could not find build directory") diff --git a/boututils/surface_average.py b/boututils/surface_average.py index d58e9a3..340d252 100644 --- a/boututils/surface_average.py +++ b/boututils/surface_average.py @@ -11,17 +11,16 @@ from boututils.calculus import deriv from boututils.int_func import int_func from .idl_tabulate import idl_tabulate -from bunch import bunchify -def surface_average(var, g, area=None): +def surface_average(var, grid, area=None): """Average a variable over a surface Parameters ---------- var : array_like 3-D or 4D variable to integrate (either [x,y,z] or [t,x,y,z]) - g : dict + grid : dict A dictionary of various grid quantities area : bool Average by flux-surface area = (B/Bp)*dl * R*dz @@ -35,8 +34,6 @@ def surface_average(var, g, area=None): s = np.ndim(var) - - if s == 4 : nx = np.shape(var)[1] ny = np.shape(var)[2] @@ -44,28 +41,21 @@ def surface_average(var, g, area=None): result = np.zeros((nx,nt)) for t in range (nt): - - result[:,t] = surface_average(var[t,:,:,:], g, area=area) + result[:,t] = surface_average(var[t,:,:,:], grid, area=area) return result elif s != 3 : - print("ERROR: surface_average var must be 3 or 4D") - return 0 + raise RuntimeError("ERROR: surface_average var must be 3D or 4D") - - # 3D [x,y,z] + # 3D [x,y,z] nx = np.shape(var)[0] ny = np.shape(var)[1] -# nz = np.shape(var)[2] - -# Use bunch to create grid structure - grid=bunchify(g) - # Calculate poloidal angle from grid + # Calculate poloidal angle from grid theta = np.zeros((nx,ny)) - #status = gen_surface(mesh=grid) ; Start generator + #status = gen_surface(mesh=grid) ; Start generator xi = -1 yi = np.arange(0,ny,dtype=int) last = 0 @@ -76,17 +66,17 @@ def surface_average(var, g, area=None): last = 1 dtheta = 2.*np.pi / np.float(ny) - r = grid.Rxy[xi,yi] - z = grid.Zxy[xi,yi] + r = grid['Rxy'][xi,yi] + z = grid['Zxy'][xi,yi] n = np.size(r) dl = old_div(np.sqrt( deriv(r)**2 + deriv(z)**2 ), dtheta) if area: - dA = (old_div(grid.Bxy[xi,yi],grid.Bpxy[xi,yi]))*r*dl + dA = (old_div(grid['Bxy'][xi,yi],grid['Bpxy'][xi,yi]))*r*dl A = int_func(np.arange(n),dA) theta[xi,yi] = 2.*np.pi*A/A[n-1] else: - nu = dl * (grid.Btxy[xi,yi]) / ((grid.Bpxy[xi,yi]) * r ) + nu = dl * (grid['Btxy'][xi,yi]) / ((grid['Bpxy'][xi,yi]) * r ) theta[xi,yi] = int_func(np.arange(n)*dtheta,nu) theta[xi,yi] = 2.*np.pi*theta[xi,yi] / theta[xi,yi[n-1]] diff --git a/boututils/volume_integral.py b/boututils/volume_integral.py index 9fdb083..7e5d4da 100644 --- a/boututils/volume_integral.py +++ b/boututils/volume_integral.py @@ -8,17 +8,15 @@ from past.utils import old_div import numpy as np from boututils.calculus import deriv -from bunch import bunchify - -def volume_integral(var, g, xr=False): +def volume_integral(var, grid, xr=False): """Integrate a variable over a volume Parameters ---------- var : array_like Variable to integrate - g : dict + grid : dict A dictionary of various grid quantities xr : (int, int), optional Range of x indices (default: all of x) @@ -32,9 +30,6 @@ def volume_integral(var, g, xr=False): s = np.ndim(var) - grid=bunchify(g) - - if s == 4 : # 4D [t,x,y,z] - integrate for each t nx = np.shape(var)[1] @@ -84,20 +79,20 @@ def volume_integral(var, g, xr=False): if (xi >= np.min(xr)) & (xi <= np.max(xr)) : dtheta = 2.*np.pi / np.float(ny) - r = grid.Rxy[xi,yi] - z = grid.Zxy[xi,yi] + r = grid['Rxy'][xi,yi] + z = grid['Zxy'][xi,yi] n = np.size(r) dl = old_div(np.sqrt( deriv(r)**2 + deriv(z)**2 ), dtheta) # Area of flux-surface - dA = (grid.Bxy[xi,yi]/grid.Bpxy[xi,yi]*dl) * (r*2.*np.pi) + dA = (grid['Bxy'][xi,yi]/grid['Bpxy'][xi,yi]*dl) * (r*2.*np.pi) # Volume if xi == nx-1 : - dpsi = (grid.psixy[xi,yi] - grid.psixy[xi-1,yi]) + dpsi = (grid['psixy'][xi,yi] - grid['psixy'][xi-1,yi]) else: - dpsi = (grid.psixy[xi+1,yi] - grid.psixy[xi,yi]) + dpsi = (grid['psixy'][xi+1,yi] - grid['psixy'][xi,yi]) - dV = dA * dpsi / (r*(grid.Bpxy[xi,yi])) # May need factor of 2pi + dV = dA * dpsi / (r*(grid['Bpxy'][xi,yi])) # May need factor of 2pi dV = np.abs(dV) result = result + np.sum(var[xi,yi] * dV) From f599a57d38e60f342cd7c4c932a83139054dbf10 Mon Sep 17 00:00:00 2001 From: David Date: Sun, 20 Sep 2020 08:51:57 +0200 Subject: [PATCH 2/2] remove bunch --- requirements.txt | 1 - setup.py | 1 - 2 files changed, 2 deletions(-) diff --git a/requirements.txt b/requirements.txt index 1dbb35b..7f6d297 100644 --- a/requirements.txt +++ b/requirements.txt @@ -3,7 +3,6 @@ matplotlib==3.2.1 h5py==2.10.0 scipy==1.4.1 future==0.18.2 -bunch==1.0.1 PyQt5==5.12.2 [mayavi] mayavi==4.6.2 [mayavi] netCDF4==1.5.3 diff --git a/setup.py b/setup.py index 260d573..5de0663 100644 --- a/setup.py +++ b/setup.py @@ -50,7 +50,6 @@ 'h5py', 'boututils', 'future', - 'bunch', 'netCDF4'], extras_require={ 'mayavi': ['mayavi', 'PyQt5']},