diff --git a/backends/cbackend.py b/backends/cbackend.py index 00d69957fe301feed21162cc8e9f9c57f139eca7..0997ff97e44aa2fe44406e99462dc7382d0e8560 100644 --- a/backends/cbackend.py +++ b/backends/cbackend.py @@ -1,4 +1,7 @@ import sympy as sp + +from pystencils.bitoperations import xor, rightShift, leftShift + try: from sympy.utilities.codegen import CCodePrinter except ImportError: @@ -210,6 +213,12 @@ class CustomSympyPrinter(CCodePrinter): if expr.func == castFunc: arg, type = expr.args return "*((%s)(& %s))" % (PointerType(type), self._print(arg)) + elif expr.func == xor: + return "(%s ^ %s" % (self._print(expr.args[0]), self._print(expr.args[1])) + elif expr.func == rightShift: + return "(%s >> %s)" % (self._print(expr.args[0]), self._print(expr.args[1])) + elif expr.func == leftShift: + return "(%s << %s)" % (self._print(expr.args[0]), self._print(expr.args[1])) else: return super(CustomSympyPrinter, self)._print_Function(expr) diff --git a/bitoperations.py b/bitoperations.py new file mode 100644 index 0000000000000000000000000000000000000000..fa73e218f1a211324c618332972a1228c4d07d31 --- /dev/null +++ b/bitoperations.py @@ -0,0 +1,5 @@ +import sympy as sp +xor = sp.Function("⊻") +rightShift = sp.Function("rshift") +leftShift = sp.Function("lshift") + diff --git a/datahandling/datahandling_interface.py b/datahandling/datahandling_interface.py index d4d1f428c71525bceaa363e6a0ed31d5227d33a0..25914e990e048b0a0984b945205c4896173b77a4 100644 --- a/datahandling/datahandling_interface.py +++ b/datahandling/datahandling_interface.py @@ -241,12 +241,13 @@ class DataHandling(ABC): def __str__(self): result = "" - rowFormat = "{:>35}|{:>21}|{:>21}\n" - separatorLine = "" * (25 + 21 + 21 + 2) + "\n" + firstColumnWidth = max(len("Name"), max(len(a) for a in self.arrayNames)) + rowFormat = "{:>%d}|{:>21}|{:>21}\n" % (firstColumnWidth,) + separatorLine = "-" * (firstColumnWidth + 21 + 21 + 2) + "\n" result += rowFormat.format("Name", "Inner (min/max)", "WithGl (min/max)") result += separatorLine for arrName in sorted(self.arrayNames): - innerMinMax = (self.min(arrName, ghostLayers=False), self.max(arrName, ghostLayers=False)) + innerMinMax = (self.min(arrName, ghostLayers=False), self.max(arrName, ghostLayers=False)) withGlMinMax = (self.min(arrName, ghostLayers=True), self.max(arrName, ghostLayers=True)) innerMinMax = "({0[0]:3.3g},{0[1]:3.3g})".format(innerMinMax) withGlMinMax = "({0[0]:3.3g},{0[1]:3.3g})".format(withGlMinMax) diff --git a/datahandling/parallel_datahandling.py b/datahandling/parallel_datahandling.py index 2b1007434a22dfdb36ff5fc20dd156b4b0011385..2113ed6b458c0e80d0039f4e2b3c9b3954afb19d 100644 --- a/datahandling/parallel_datahandling.py +++ b/datahandling/parallel_datahandling.py @@ -84,8 +84,6 @@ class ParallelDataHandling(DataHandling): ghostLayers = self.defaultGhostLayers if layout is None: layout = self.defaultLayout - if latexName is None: - latexName = name if len(self.blocks) == 0: raise ValueError("Data handling expects that each process has at least one block") if hasattr(dtype, 'type'): @@ -118,13 +116,14 @@ class ParallelDataHandling(DataHandling): if indexDimensions == 1: shape += (fSize, ) - assert all(f.name != latexName for f in self.fields.values()), "Symbolic field with this name already exists" + assert all(f.name != name for f in self.fields.values()), "Symbolic field with this name already exists" - self.fields[name] = Field.createGeneric(latexName, self.dim, dtype, indexDimensions, layout, + self.fields[name] = Field.createGeneric(name, self.dim, dtype, indexDimensions, layout, indexShape=(fSize,) if indexDimensions > 0 else None) - self._fieldNameToCpuDataName[latexName] = name + self.fields[name].latexName = latexName + self._fieldNameToCpuDataName[name] = name if gpu: - self._fieldNameToGpuDataName[latexName] = self.GPU_DATA_PREFIX + name + self._fieldNameToGpuDataName[name] = self.GPU_DATA_PREFIX + name return self.fields[name] def hasData(self, name): @@ -139,7 +138,7 @@ class ParallelDataHandling(DataHandling): return tuple(self._customDataNames) def addArrayLike(self, name, nameOfTemplateField, latexName=None, cpu=True, gpu=False): - self.addArray(name, latexName=latexName, cpu=cpu, gpu=gpu, **self._fieldInformation[nameOfTemplateField]) + return self.addArray(name, latexName=latexName, cpu=cpu, gpu=gpu, **self._fieldInformation[nameOfTemplateField]) def swap(self, name1, name2, gpu=False): if gpu: @@ -178,16 +177,21 @@ class ParallelDataHandling(DataHandling): if sliceObj is None: sliceObj = tuple([slice(None, None, None)] * self.dim) if self.dim == 2: - sliceObj += (0.5,) + sliceObj = sliceObj[:2] + (0.5,) + sliceObj[2:] - array = wlb.field.gatherField(self.blocks, name, sliceObj, allGather) + lastElement = sliceObj[3:] + + array = wlb.field.gatherField(self.blocks, name, sliceObj[:3], allGather) if array is None: return None - if self.fields[name].indexDimensions == 0: - array = array[..., 0] if self.dim == 2: array = array[:, :, 0] + if lastElement and self.fields[name].indexDimensions > 0: + array = array[..., lastElement[0]] + if self.fields[name].indexDimensions == 0: + array = array[..., 0] + return array def _normalizeArrShape(self, arr, indexDimensions): diff --git a/datahandling/serial_datahandling.py b/datahandling/serial_datahandling.py index a70b48de17992d55b8fdb00ea2da252620742ec1..85970d3ad6e490d17375fd4d40310dd22b4cbd97 100644 --- a/datahandling/serial_datahandling.py +++ b/datahandling/serial_datahandling.py @@ -31,7 +31,6 @@ class SerialDataHandling(DataHandling): self.defaultGhostLayers = defaultGhostLayers self.defaultLayout = defaultLayout self._fields = DotDict() - self._fieldLatexNameToDataName = {} self.cpuArrays = DotDict() self.gpuArrays = DotDict() self.customDataCpu = DotDict() @@ -74,8 +73,6 @@ class SerialDataHandling(DataHandling): ghostLayers = self.defaultGhostLayers if layout is None: layout = self.defaultLayout - if latexName is None: - latexName = name kwargs = { 'shape': tuple(s + 2 * ghostLayers for s in self._domainSize), @@ -109,10 +106,10 @@ class SerialDataHandling(DataHandling): raise ValueError("GPU Field with this name already exists") self.gpuArrays[name] = gpuarray.to_gpu(cpuArr) - assert all(f.name != latexName for f in self.fields.values()), "Symbolic field with this name already exists" - self.fields[name] = Field.createFixedSize(latexName, shape=kwargs['shape'], indexDimensions=indexDimensions, + assert all(f.name != name for f in self.fields.values()), "Symbolic field with this name already exists" + self.fields[name] = Field.createFixedSize(name, shape=kwargs['shape'], indexDimensions=indexDimensions, dtype=kwargs['dtype'], layout=layoutTuple) - self._fieldLatexNameToDataName[latexName] = name + self.fields[name].latexName = latexName return self.fields[name] def addCustomData(self, name, cpuCreationFunction, @@ -136,7 +133,7 @@ class SerialDataHandling(DataHandling): return name in self.fields def addArrayLike(self, name, nameOfTemplateField, latexName=None, cpu=True, gpu=False): - self.addArray(name, latexName=latexName, cpu=cpu, gpu=gpu, **self._fieldInformation[nameOfTemplateField]) + return self.addArray(name, latexName=latexName, cpu=cpu, gpu=gpu, **self._fieldInformation[nameOfTemplateField]) def iterate(self, sliceObj=None, gpu=False, ghostLayers=True, innerGhostLayers=True): if ghostLayers is True: @@ -172,12 +169,15 @@ class SerialDataHandling(DataHandling): glToRemove = 0 arr = self.cpuArrays[name] indDimensions = self.fields[name].indexDimensions + spatialDimensions = self.fields[name].spatialDimensions + arr = removeGhostLayers(arr, indexDimensions=indDimensions, ghostLayers=glToRemove) if sliceObj is not None: - sliceObj = normalizeSlice(sliceObj, arr.shape[:-indDimensions] if indDimensions > 0 else arr.shape) - sliceObj = tuple(s if type(s) is slice else slice(s, s + 1, None) for s in sliceObj) - arr = arr[sliceObj] + normalizedSlice = normalizeSlice(sliceObj[:spatialDimensions], arr.shape[:spatialDimensions]) + normalizedSlice = tuple(s if type(s) is slice else slice(s, s + 1, None) for s in normalizedSlice) + normalizedSlice += sliceObj[spatialDimensions:] + arr = arr[normalizedSlice] else: arr = arr.view() arr.flags.writeable = False @@ -198,7 +198,7 @@ class SerialDataHandling(DataHandling): self.toGpu(name) def runKernel(self, kernelFunc, *args, **kwargs): - dataUsedInKernel = [self._fieldLatexNameToDataName[p.fieldName] + dataUsedInKernel = [p.fieldName for p in kernelFunc.parameters if p.isFieldPtrArgument and p.fieldName not in kwargs] arrays = self.gpuArrays if kernelFunc.ast.backend == 'gpucuda' else self.cpuArrays arrayParams = {name: arrays[name] for name in dataUsedInKernel} @@ -270,11 +270,11 @@ class SerialDataHandling(DataHandling): if target == 'cpu': def resultFunctor(): - for func in resultFunctors: + for name, func in zip(names, resultFunctors): func(pdfs=self.cpuArrays[name]) else: def resultFunctor(): - for func in resultFunctors: + for name, func in zip(names, resultFunctors): func(pdfs=self.gpuArrays[name]) return resultFunctor diff --git a/field.py b/field.py index def8e0dd71b57988bd2b0da7e40c6c9213eb4210..b64e567e2ae07f09faf1f574fceadabc5a66bcad 100644 --- a/field.py +++ b/field.py @@ -188,6 +188,7 @@ class Field(object): self._layout = normalizeLayout(layout) self.shape = shape self.strides = strides + self.latexName = None def newFieldWithDifferentName(self, newName): return Field(newName, self.fieldType, self._dtype, self._layout, self.shape, self.strides) @@ -381,10 +382,11 @@ class Field(object): return self._offsetName def _latex(self, arg): + n = self._field.latexName if self._field.latexName else self._field.name if self._superscript: - return "{{%s}_{%s}^{%s}}" % (self._field.name, self._offsetName, self._superscript) + return "{{%s}_{%s}^{%s}}" % (n, self._offsetName, self._superscript) else: - return "{{%s}_{%s}}" % (self._field.name, self._offsetName) + return "{{%s}_{%s}}" % (n, self._offsetName) @property def index(self): diff --git a/finitedifferences.py b/finitedifferences.py index 835934bb07e278093346e40bc2ad006c8a064083..9f1152d151e9806b72497ec16d81221536829cda 100644 --- a/finitedifferences.py +++ b/finitedifferences.py @@ -247,7 +247,7 @@ class Transient(sp.Function): @property def scalar(self): if self.scalarIndex is None: - return self.args[0].field + return self.args[0].field.center else: return self.args[0].field(self.scalarIndex) @@ -310,9 +310,9 @@ class Discretization2ndOrder: if len(solveResult) != 1: raise ValueError("Could not solve for transient term" + str(solveResult)) rhs = solveResult.pop() - idx = 0 if transientTerm.scalarIndex is None else transientTerm.scalarIndex # explicit euler - return transientTerm.scalar(idx) + self.dt * self._discretizeSpatial(rhs) + return transientTerm.scalar + self.dt * self._discretizeSpatial(rhs) else: + print(transientTerms) raise NotImplementedError("Cannot discretize expression with more than one transient term")