diff --git a/.arclint b/.arclint
new file mode 100644
index 0000000..73c2265
--- /dev/null
+++ b/.arclint
@@ -0,0 +1,18 @@
+{
+ "linters": {
+ "python": {
+ "type": "flake8",
+ "include": [ "(\\.py$)", "(^hbm$)" ],
+ "exclude": ["(^third-party/)", "(legacy/)"],
+ "severity.rules": {
+ "(^E)": "error",
+ "(^W)": "warning",
+ "(^F)": "advice"
+ }
+ },
+
+ "merges": {
+ "type": "merge-conflict"
+ }
+ }
+}
diff --git a/.flake8 b/.flake8
new file mode 100644
index 0000000..15fc7e3
--- /dev/null
+++ b/.flake8
@@ -0,0 +1,2 @@
+[flake8]
+max-line-length = 80
diff --git a/doc/Makefile b/doc/Makefile
index 0c5906e..ed502aa 100644
--- a/doc/Makefile
+++ b/doc/Makefile
@@ -1,32 +1,32 @@
DOXYGEN = doxygen
SPHINX = sphinx-build
SPHINX_SOURCE_DIR = sphinx/source
SPHINX_SOURCES := $(wildcard $(SPHINX_SOURCE_DIR)/*.rst)
SPHINX_SOURCES += $(SPHINX_SOURCE_DIR)/conf.py
BUILD_DIR = build
DOXYGEN_BUILD_DIR = $(BUILD_DIR)/doxygen
SPHINX_BUILD_DIR = $(BUILD_DIR)/sphinx/html
# Default target
-all:sphinx
+all:doxygen sphinx
# Phony targets for aliases
.PHONY:doxygen sphinx clean
# Target aliases
doxygen:$(DOXYGEN_BUILD_DIR)
sphinx:$(SPHINX_BUILD_DIR)
# Build doxygen documentation
$(DOXYGEN_BUILD_DIR):doxygen/Doxyfile
@mkdir -p $@
$(DOXYGEN) $<
$(SPHINX_BUILD_DIR):$(SPHINX_SOURCES) $(DOXYGEN_BUILD_DIR)
@mkdir -p $@
$(SPHINX) -b html $(SPHINX_SOURCE_DIR) $@
clean:
$(RM) -r $(BUILD_DIR)
diff --git a/examples/SConstruct b/examples/SConstruct
deleted file mode 100644
index af60e3b..0000000
--- a/examples/SConstruct
+++ /dev/null
@@ -1,31 +0,0 @@
-# -*- mode:python; coding: utf-8 -*-
-# vim: set ft=python:
-
-# @file
-# @section LICENSE
-#
-# Copyright (©) 2016-19 EPFL (École Polytechnique Fédérale de Lausanne),
-# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
-#
-# This program is free software: you can redistribute it and/or modify
-# it under the terms of the GNU Affero General Public License as published
-# by the Free Software Foundation, either version 3 of the License, or
-# (at your option) any later version.
-#
-# This program is distributed in the hope that it will be useful,
-# but WITHOUT ANY WARRANTY; without even the implied warranty of
-# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
-# GNU Affero General Public License for more details.
-#
-# You should have received a copy of the GNU Affero General Public License
-# along with this program. If not, see .
-
-import os
-
-env = Environment(CC = 'g++',
- CXXFLAGS = ['-std=c++11', '-g', '-fopenmp', '-O3'],
- TOOLS = ['default', 'tamaas'],
- TAMAAS_PATH = '../build-release/',
- ENV = os.environ)
-
-env.Program('rough_surface.cc')
diff --git a/examples/rough_contact.py b/examples/rough_contact.py
index 10f145b..abfb80e 100644
--- a/examples/rough_contact.py
+++ b/examples/rough_contact.py
@@ -1,74 +1,74 @@
#!/usr/bin/env python3
# @file
# @section LICENSE
#
# Copyright (©) 2016-19 EPFL (École Polytechnique Fédérale de Lausanne),
# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU Affero General Public License as published
# by the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Affero General Public License for more details.
#
# You should have received a copy of the GNU Affero General Public License
# along with this program. If not, see .
import tamaas as tm
import matplotlib.pyplot as plt
# Initialize threads and fftw
tm.initialize()
tm.set_log_level(tm.LogLevel.info) # Show progression of solver
# Surface size
n = 512
# Surface generator
sg = tm.SurfaceGeneratorFilter2D()
sg.setSizes([n, n])
sg.random_seed = 0
# Spectrum
spectrum = tm.Isopowerlaw2D()
sg.setFilter(spectrum)
# Parameters
spectrum.q0 = 16
spectrum.q1 = 16
spectrum.q2 = 64
spectrum.hurst = 0.8
# Generating surface
surface = sg.buildSurface()
-#surface /= tm.SurfaceStatistics.computeSpectralRMSSlope(surface)
+# surface /= tm.SurfaceStatistics.computeSpectralRMSSlope(surface)
surface /= n
-#print(spectrum.rmsSlopes())
-#print(tm.SurfaceStatistics.computeRMSSlope(surface))
+# print(spectrum.rmsSlopes())
+# print(tm.SurfaceStatistics.computeRMSSlope(surface))
plt.imshow(surface)
# Creating model
model = tm.ModelFactory.createModel(tm.model_type_basic_2d, [1., 1.], [n, n])
# Solver
solver = tm.PolonskyKeerRey(model, surface, 1e-12,
tm.PolonskyKeerRey.pressure,
tm.PolonskyKeerRey.pressure)
# Solve for target pressure
p_target = 0.1
solver.solve(p_target)
plt.figure()
plt.imshow(model.getTraction())
print(model.getTraction().mean())
# Cleanup threads
tm.finalize()
plt.show()
diff --git a/examples/site_scons/site_init.py b/examples/site_scons/site_init.py
deleted file mode 100644
index a6745c5..0000000
--- a/examples/site_scons/site_init.py
+++ /dev/null
@@ -1,151 +0,0 @@
-# @file
-# @section LICENSE
-#
-# Copyright (©) 2016-19 EPFL (École Polytechnique Fédérale de Lausanne),
-# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
-#
-# This program is free software: you can redistribute it and/or modify
-# it under the terms of the GNU Affero General Public License as published
-# by the Free Software Foundation, either version 3 of the License, or
-# (at your option) any later version.
-#
-# This program is distributed in the hope that it will be useful,
-# but WITHOUT ANY WARRANTY; without even the implied warranty of
-# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
-# GNU Affero General Public License for more details.
-#
-# You should have received a copy of the GNU Affero General Public License
-# along with this program. If not, see .
-
-from SCons.Script import *
-from os.path import *
-
-import os
-
-
-def fftw(env):
- """A Tool to search for fftw headers and libraries"""
-
- if env.GetOption('clean'):
- return
- env.SetDefault(FFTW_VERSION='3')
- env.SetDefault(FFTW_LIBRARY_WISH=[])
-
- print("Building for fftw version {}".format(env['FFTW_VERSION']))
-
- if 'FFTW_INCLUDE_DIR' in env:
- env.AppendUnique(CPPPATH=env['FFTW_INCLUDE_DIR'])
-
- if 'FFTW_LIBRARY_DIR' in env:
- env.AppendUnique(LIBPATH=env['FFTW_LIBRARY_DIR'])
-
- version = env['FFTW_VERSION']
- if version == "2":
- lib_names = {'main': 'fftw'}
- inc_names = ['fftw.h']
- else:
- lib_names = {'main': 'fftw3',
- 'thread': 'fftw3_threads',
- 'omp': 'fftw3_omp'}
- inc_names = ['fftw3.h']
-
- try:
- libs = [lib_names[i] for i in env['FFTW_LIBRARY_WISH']]
- except:
- raise SCons.Errors.StopError(
- 'Incompatible wishlist {0} from version {1}'.format(
- env['FFTW_LIBRARY_WISH'], env['FFTW_VERSION']))
-
- env.AppendUnique(LIBS=libs)
- if version == "2":
- env.Append(LIBS='m')
-
- conf = Configure(env)
- if not conf.CheckLibWithHeader(libs, inc_names, 'c++'):
- raise SCons.Errors.StopError(
- 'Failed to find libraries {0} or '
- 'headers {1}.'.format(str(lib_names), str(inc_names)))
-
- env = conf.Finish()
-
-
-def criterion(env):
- """A tool to configure for Criterion (test suite)"""
- if env.GetOption('clean'):
- return
-
- env.AppendUnique(CPPPATH=[Dir("#third-party/Criterion/include")],
- LIBPATH=[abspath(str(Dir("#third-party/Criterion/build")))])
- conf = Configure(env)
- # import pdb
- # pdb.set_trace()
- if not conf.CheckLibWithHeader('criterion',
- 'criterion/criterion.h',
- 'c++'):
- raise SCons.Errors.StopError(
- 'Failed to find library Criterion or criterion.h header in third-party')
- env = conf.Finish()
-
-
-def boost(env):
- """A tool to confugure for boost"""
- if env.GetOption('clean'):
- return
-
- if 'BOOST_INCLUDE_DIR' in env:
- env.AppendUnique(CPPPATH=env['BOOST_INCLUDE_DIR'])
- conf = Configure(env)
- if not conf.CheckCXXHeader('boost/version.hpp'):
- raise SCons.Errors.StopError(
- 'Failed to find Boost library')
- env = conf.Finish()
-
-
-def thrust(env):
- """A tool to confugure for thrust"""
- if env.GetOption('clean'):
- return
-
- if 'THRUST_INCLUDE_DIR' in env:
- env.AppendUnique(CPPPATH=env['THRUST_INCLUDE_DIR'])
- conf = Configure(env)
- if not conf.CheckCXXHeader('thrust/version.h'):
- raise SCons.Errors.StopError(
- 'Failed to find Thrust library')
- env = conf.Finish()
-
-
-def tamaas(env):
- """Sets correct environement variables for Tamaas.
- - 'TAMAAS_PATH' is list of potential tamaas repositories
- - 'TAMAAS_BUILD_TYPE' is the build type of tamaas (release, debug, profiling)"""
- if env.GetOption("clean"):
- return
-
- try:
- build_type = env['TAMAAS_BUILD_TYPE']
- except:
- build_type = "release"
-
- if 'TAMAAS_PATH' not in env:
- print("Please set the 'TAMAAS_PATH' variable in your scons environment")
- Exit(1)
-
- tm_path = env['TAMAAS_PATH']
- include_dir = tm_path + '/src'
- subdirs = "core model solvers gpu bem surface".split(" ")
- include_dirs = [abspath(join(include_dir, x)) for x in subdirs]
- lib_dir = join(tm_path, 'build-{}/src'.format(build_type))
-
- env.AppendUnique(CPPPATH = include_dirs)
- env.AppendUnique(CPPPATH = ['/opt/cuda/include']) # TODO fix this dirty hack
- env.AppendUnique(LIBPATH = [abspath(lib_dir)])
- env.AppendUnique(RPATH = [abspath(lib_dir)])
- env.AppendUnique(CXXFLAGS = ["-std=c++14"])
-
- conf = Configure(env)
- if not conf.CheckLibWithHeader("Tamaas", "tamaas.hh", language="CXX"):
- raise SCons.Errors.StopError(
- "Tamaas ({}) not found in {}".format(build_type, abspath(tm_path)))
- env['TAMAAS_PATH'] = abspath(tm_path)
- env = conf.Finish()
diff --git a/examples/statistics.py b/examples/statistics.py
index 1168807..528c00a 100644
--- a/examples/statistics.py
+++ b/examples/statistics.py
@@ -1,86 +1,87 @@
#!/usr/bin/env python3
# @file
# @section LICENSE
#
# Copyright (©) 2016-19 EPFL (École Polytechnique Fédérale de Lausanne),
# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU Affero General Public License as published
# by the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Affero General Public License for more details.
#
# You should have received a copy of the GNU Affero General Public License
# along with this program. If not, see .
import os
import numpy as np
import matplotlib.pyplot as plt
import tamaas as tm
from matplotlib.colors import LogNorm
# Initialize threads and fftw
tm.initialize()
tm.set_log_level(tm.LogLevel.info) # Show progression of solver
# Surface size
n = 512
# Surface generator
sg = tm.SurfaceGeneratorFilter2D()
sg.setSizes([n, n])
sg.random_seed = 0
# Spectrum
spectrum = tm.Isopowerlaw2D()
sg.setFilter(spectrum)
# Parameters
spectrum.q0 = 16
spectrum.q1 = 16
spectrum.q2 = 64
spectrum.hurst = 0.8
# Generating surface
surface = sg.buildSurface()
# Computing PSD and shifting for plot
psd = tm.Statistics2D.computePowerSpectrum(surface)
psd = np.fft.fftshift(psd, axes=0)
plt.imshow(psd.real, norm=LogNorm())
plt.gca().set_title('Power Spectrum Density')
plt.gcf().tight_layout()
# Computing autocorrelation and shifting for plot
acf = tm.Statistics2D.computeAutocorrelation(surface)
acf = np.fft.fftshift(acf)
plt.figure()
plt.imshow(acf)
plt.gca().set_title('Autocorrelation')
plt.gcf().tight_layout()
plt.show()
# Write the rough surface in paraview's VTK format
try:
import uvw
try:
os.mkdir('paraview')
- except FileExistsError as e:
+ except FileExistsError:
pass
x = np.linspace(0, 1, n, endpoint=True)
y = x.copy()
- with uvw.RectilinearGrid(os.path.join('paraview', 'surface.vtr'), (x, y)) as grid:
+ with uvw.RectilinearGrid(os.path.join('paraview', 'surface.vtr'),
+ (x, y)) as grid:
grid.addPointData(uvw.DataArray(surface, range(2), 'surface'))
-except ImportError as e:
+except ImportError:
print("uvw not installed, will not write VTK file")
pass
diff --git a/python/tamaas/__init__.py b/python/tamaas/__init__.py
index 00e277c..7818ba5 100644
--- a/python/tamaas/__init__.py
+++ b/python/tamaas/__init__.py
@@ -1,71 +1,73 @@
# @file
# @section LICENSE
#
# Copyright (©) 2016-19 EPFL (École Polytechnique Fédérale de Lausanne),
# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU Affero General Public License as published
# by the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Affero General Public License for more details.
#
# You should have received a copy of the GNU Affero General Public License
# along with this program. If not, see .
"""
Tamaas is a library dedicated to the fast treatment of rough contact problems.
See README.md and examples for guidance as how to use Tamaas.
See __authors__, __license__, __copyright__ for extra information about Tamaas.
"""
__authors__ = [
"Guillaume Anciaux ",
"Lucas Frérot ",
"Son Pham-Ba ",
"Valentine Rey ",
]
__copyright__ = "Copyright (©) 2016-19 EPFL " \
+ "(École Polytechnique Fédérale de Lausanne), " \
+ "Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)"
__license__ = "AGPLv3"
__maintainer__ = "Lucas Frérot"
__email__ = "lucas.frerot@epfl.ch"
-try:
- from ._tamaas import * # noqa
+try:
+ from . import _tamaas as _tm
# Redefinition of model_type constants (for compatibility)
- model_type_basic_1d = model_type.basic_1d
- model_type_basic_2d = model_type.basic_2d
- model_type_surface_1d = model_type.surface_1d
- model_type_surface_2d = model_type.surface_2d
- model_type_volume_1d = model_type.volume_1d
- model_type_volume_2d = model_type.volume_2d
+ model_type_basic_1d = _tm.model_type.basic_1d
+ model_type_basic_2d = _tm.model_type.basic_2d
+ model_type_surface_1d = _tm.model_type.surface_1d
+ model_type_surface_2d = _tm.model_type.surface_2d
+ model_type_volume_1d = _tm.model_type.volume_1d
+ model_type_volume_2d = _tm.model_type.volume_2d
+
+ from ._tamaas import * # noqa
except ImportError as e:
print("Error trying to import _tamaas:\n{}".format(e))
raise e
class ParamHelper:
"""Legacy class to manage parameters/setters/getters"""
def __init__(self, obj):
self.obj = obj
def set(self, params):
for key in params:
setter_name = 'set' + key
try:
accessor = getattr(self.obj, setter_name)
accessor(params[key])
except Exception:
print("Setter '{}({})' does not exist for object {}"
.format(setter_name, type(params[key]), self.obj))
diff --git a/site_scons/site_tools/nvcc.py b/site_scons/site_tools/nvcc.py
index a385268..7751922 100644
--- a/site_scons/site_tools/nvcc.py
+++ b/site_scons/site_tools/nvcc.py
@@ -1,205 +1,223 @@
# -*- coding: utf-8 -*-
# @file
# @section LICENSE
#
# Copyright (©) 2016-19 EPFL (École Polytechnique Fédérale de Lausanne),
# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU Affero General Public License as published
# by the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Affero General Public License for more details.
#
# You should have received a copy of the GNU Affero General Public License
# along with this program. If not, see .
"""
SCons tool to setup compilation environment for NVidia's NVCC
Based on version
https://github.com/thrust/thrust/blob/master/site_scons/site_tools/nvcc.py
"""
from __future__ import print_function
import SCons
from SCons.Errors import StopError
import platform
from os.path import join
known_components = {
'cufft': 'cufft.h',
'cufftw': 'cufftw.h'
}
def _find_cuda_paths(env):
"""Finding cuda paths"""
if 'CUDA_TOOLKIT_PATH' not in env:
raise StopError("CUDA_TOOLKIT_PATH variable not found in environment")
cuda_path = env['CUDA_TOOLKIT_PATH']
inc_path = join(cuda_path, 'include')
lib_path = join(cuda_path, 'lib')
bin_path = join(cuda_path, 'bin')
if platform.machine()[-2:] == '64':
lib_path += '64'
# add nvcc's location to PATH
env.PrependENVPath('PATH', bin_path)
# add the location of the cudart shared library to LD_LIBRARY_PATH as well
# this allows us to execute CUDA programs during the build
env.PrependENVPath('LD_LIBRARY_PATH', lib_path)
# add to include paths
env.AppendUnique(CPPPATH=[inc_path])
# add to library paths
env.AppendUnique(LIBPATH=[lib_path])
def _configure_cuda_components(env):
"""Finding required components in cuda"""
if env.GetOption('clean'):
return
if 'CUDA_COMPONENTS' in env:
for component in env['CUDA_COMPONENTS']:
if component not in known_components:
raise StopError("Unknown cuda component '{}'".format(component))
conf = SCons.Script.Configure(env)
- if not conf.CheckLibWithHeader(component, known_components[component], 'c++'):
+ if not conf.CheckLibWithHeader(component,
+ known_components[component],
+ 'c++'):
raise StopError("Failed to find library {} or header {}".format(
component, known_components[component]))
env = conf.Finish()
def _define_compile_cuda(env):
"""Define the compile string for files that need nvcc"""
env['NVCC'] = env.Detect('nvcc')
# Setting compilation command
env['NVCCCOM'] = \
- '$NVCC $CUDA_ARCH_FLAG -o $TARGET -c $NVCC_CXXFLAGS $NVCC_CCFLAGS $_CCCOMCOM -x cu $SOURCES'
+ '$NVCC $CUDA_ARCH_FLAG -o $TARGET -c $NVCC_CXXFLAGS $NVCC_CCFLAGS ' \
+ + '$_CCCOMCOM -x cu $SOURCES'
env['SHNVCCCOM'] = \
- '$NVCC $CUDA_ARCH_FLAG -o $TARGET -dc -Xcompiler -fPIC $NVCC_CXXFLAGS $NVCC_CCFLAGS $_CCCOMCOM -x cu $SOURCES'
+ '$NVCC $CUDA_ARCH_FLAG -o $TARGET -dc -Xcompiler -fPIC ' \
+ + '$NVCC_CXXFLAGS $NVCC_CCFLAGS $_CCCOMCOM -x cu $SOURCES'
# Constructing proper cc flags
env['_NVCC_BARE_CCFLAGS'] = \
'${_concat("", CCFLAGS, "", __env__, _NVCC_BARE_FILTER)}'
env['_NVCC_PASSED_CCFLAGS'] = \
- '${_concat("-Xcompiler ", CCFLAGS, "", __env__, _NVCC_COMPILER_PASSED_FILTER)}'
+ '${_concat("-Xcompiler ", CCFLAGS, "", __env__, ' \
+ + '_NVCC_COMPILER_PASSED_FILTER)}'
# Constructing proper cxx flags
env['_NVCC_BARE_CXXFLAGS'] = \
'${_concat("", CXXFLAGS, "", __env__, _NVCC_BARE_FILTER)}'
env['_NVCC_PASSED_CXXFLAGS'] = \
- '${_concat("-Xcompiler ", CXXFLAGS, "", __env__, _NVCC_COMPILER_PASSED_FILTER)}'
+ '${_concat("-Xcompiler ", CXXFLAGS, "", __env__, ' \
+ + '_NVCC_COMPILER_PASSED_FILTER)}'
# Putting all together
env['NVCC_CCFLAGS'] = '$_NVCC_BARE_CCFLAGS $_NVCC_PASSED_CCFLAGS'
env['NVCC_CXXFLAGS'] = '$_NVCC_BARE_CXXFLAGS $_NVCC_PASSED_CXXFLAGS'
def _define_link_cuda(env):
"""Define the link string"""
# Fixing rpaths
env['RPATHPREFIX'] = '-rpath='
env['__RPATH'] = '${_concat("-Xlinker ", _RPATH, "", __env__)}'
env['LINK'] = env['NVCC']
env['SHLINK'] = env['NVCC']
# Replacing old link command strings
env['LINKCOM'] = \
- '$LINK $CUDA_ARCH_FLAG -o $TARGET $NVCC_LINKFLAGS $__RPATH $SOURCES $_LIBDIRFLAGS $_LIBFLAGS'
+ '$LINK $CUDA_ARCH_FLAG -o $TARGET $NVCC_LINKFLAGS $__RPATH ' \
+ + '$SOURCES $_LIBDIRFLAGS $_LIBFLAGS'
env['SHLINKCOM'] = \
- '$SHLINK -shared $CUDA_ARCH_FLAG -o $TARGET $NVCC_SHLINKFLAGS $__SHLIBVERSIONFLAGS $__RPATH $SOURCES $_LIBDIRFLAGS $_LIBFLAGS'
+ '$SHLINK -shared $CUDA_ARCH_FLAG -o $TARGET $NVCC_SHLINKFLAGS ' \
+ + '$__SHLIBVERSIONFLAGS $__RPATH $SOURCES $_LIBDIRFLAGS $_LIBFLAGS'
# Constructing proper static linker flags
env['_NVCC_BARE_LINKFLAGS'] = \
'${_concat("", LINKFLAGS, "", __env__, _NVCC_BARE_FILTER)}'
env['_NVCC_COMPILER_LINKFLAGS'] = \
- '${_concat("-Xcompiler ", LINKFLAGS, "", __env__, _NVCC_COMPILER_FILTER)}'
+ '${_concat("-Xcompiler ", LINKFLAGS, "", __env__, ' \
+ + '_NVCC_COMPILER_FILTER)}'
env['_NVCC_PASSED_LINKFLAGS'] = \
- '${_concat("-Xlinker ", LINKFLAGS, "", __env__, _NVCC_LINK_PASSED_FILTER)}'
+ '${_concat("-Xlinker ", LINKFLAGS, "", __env__, ' \
+ + '_NVCC_LINK_PASSED_FILTER)}'
# Constructing proper shared linker flags
env['_NVCC_BARE_SHLINKFLAGS'] = \
'${_concat("", SHLINKFLAGS, "", __env__, _NVCC_BARE_FILTER)}'
env['_NVCC_COMPILER_SHLINKFLAGS'] = \
- '${_concat("-Xcompiler ", LINKFLAGS, "", __env__, _NVCC_COMPILER_FILTER)}'
+ '${_concat("-Xcompiler ", LINKFLAGS, "", __env__, ' \
+ + '_NVCC_COMPILER_FILTER)}'
env['_NVCC_PASSED_SHLINKFLAGS'] = \
- '${_concat("-Xlinker ", SHLINKFLAGS, "", __env__, _NVCC_LINK_PASSED_FILTER)}'
+ '${_concat("-Xlinker ", SHLINKFLAGS, "", __env__, ' \
+ + '_NVCC_LINK_PASSED_FILTER)}'
# Putting all together
env['NVCC_LINKFLAGS'] = \
- '$_NVCC_BARE_LINKFLAGS $_NVCC_COMPILER_LINKFLAGS $_NVCC_PASSED_LINKFLAGS'
+ '$_NVCC_BARE_LINKFLAGS $_NVCC_COMPILER_LINKFLAGS ' \
+ + '$_NVCC_PASSED_LINKFLAGS'
env['NVCC_SHLINKFLAGS'] = env['NVCC_LINKFLAGS']
env['NVCC_SHLINKFLAGS'] += \
- ' $_NVCC_BARE_SHLINKFLAGS $_NVCC_COMPILER_SHLINKFLAGS $_NVCC_PASSED_SHLINKFLAGS'
+ ' $_NVCC_BARE_SHLINKFLAGS $_NVCC_COMPILER_SHLINKFLAGS ' \
+ + '$_NVCC_PASSED_SHLINKFLAGS'
def _define_commands(env):
"""Defining the command strings"""
# Flags allowed by nvcc
bare_flags = """-std=c++11 -O0 -O1 -O2 -O3 -g -pg -G -w""".split()
# Experimental flags
bare_flags += """-expt-extended-lambda -expt-relaxed-constexpr""".split()
# Flags that must be passed to compiler at link time
compiler_flags = """-fopenmp -qopenmp""".split()
# Pass flags bare to nvcc
- env['_NVCC_BARE_FILTER'] = lambda flags: list(filter(lambda flag: flag in bare_flags, flags))
+ env['_NVCC_BARE_FILTER'] = lambda flags: list(filter(
+ lambda flag: flag in bare_flags, flags))
# Must prepend -Xcompiler, even at link time
env['_NVCC_COMPILER_FILTER'] = lambda flags: list(filter(
lambda flag: flag in compiler_flags, flags))
# Prepend -Xcompiler
env['_NVCC_COMPILER_PASSED_FILTER'] = lambda flags: list(filter(
lambda flag: flag not in set(bare_flags), flags))
# Prepend -Xlinker
env['_NVCC_LINKED_PASSED_FILTER'] = lambda flags: list(filter(
lambda flag: flag not in set(bare_flags) | set(compiler_flags), flags))
_define_compile_cuda(env)
_define_link_cuda(env)
def _add_actions_cuda(env):
- """Adding actions to .cu files to compile with nvcc. Other files are not affected."""
+ """
+ Adding actions to .cu files to compile with nvcc.
+ Other files are not affected.
+ """
nvcc_action = SCons.Action.Action('$NVCCCOM', '$NVCCCOMSTR')
shnvcc_action = SCons.Action.Action('$SHNVCCCOM', '$NVCCCOMSTR')
static, shared = SCons.Tool.createObjBuilders(env)
# Compiling with nvcc action added to detected .cu files
static.add_action('.cu', nvcc_action)
shared.add_action('.cu', shnvcc_action)
# Emitter to qualify correctly object code type
static.add_emitter('.cu', SCons.Defaults.StaticObjectEmitter)
shared.add_emitter('.cu', SCons.Defaults.SharedObjectEmitter)
env['CXXCOM'] = '$NVCCCOM'
env['SHCXXCOM'] = '$SHNVCCCOM'
# Scanner for dependency calculations
SCons.Tool.SourceFileScanner.add_scanner('.cu', SCons.Scanner.C.CScanner())
def generate(env):
"""Setup environment for compiling .cu files with nvcc"""
_find_cuda_paths(env)
_add_actions_cuda(env)
_define_commands(env)
_configure_cuda_components(env)
def exists(env):
return env.Detect('nvcc')
diff --git a/site_scons/site_tools/tamaas.py b/site_scons/site_tools/tamaas.py
index 6c6065e..d018941 100644
--- a/site_scons/site_tools/tamaas.py
+++ b/site_scons/site_tools/tamaas.py
@@ -1,61 +1,63 @@
# -*- coding: utf-8 -*-
# @file
# @section LICENSE
#
# Copyright (©) 2016-19 EPFL (École Polytechnique Fédérale de Lausanne),
# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU Affero General Public License as published
# by the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Affero General Public License for more details.
#
# You should have received a copy of the GNU Affero General Public License
# along with this program. If not, see .
import os
from SCons.Errors import StopError
from SCons.Script import Configure
def generate(env):
"""Sets correct environement variables for Tamaas.
- 'TAMAAS_PATH' is the path to the build directory of Tamaas"""
if env.GetOption("clean"):
return
if 'TAMAAS_PATH' not in env:
raise StopError("Please set the 'TAMAAS_PATH' variable in your " +
"scons environment")
tm_path = env['TAMAAS_PATH']
include_dir = tm_path + '/src'
subdirs = "core model solvers bem surface percolation".split(" ")
- include_dirs = [os.path.abspath(os.path.join(include_dir, x)) for x in subdirs]
+ include_dirs = [os.path.abspath(os.path.join(include_dir, x))
+ for x in subdirs]
lib_dir = include_dir
env.AppendUnique(CPPPATH=include_dirs)
env.AppendUnique(LIBPATH=[os.path.abspath(lib_dir)])
env.AppendUnique(RPATH=[os.path.abspath(lib_dir)])
env.AppendUnique(CXXFLAGS=["-std=c++14"])
- env.AppendUnique(CPPDEFINES=["THRUST_DEVICE_SYSTEM=THRUST_DEVICE_SYSTEM_OMP"])
+ env.AppendUnique(
+ CPPDEFINES=["THRUST_DEVICE_SYSTEM=THRUST_DEVICE_SYSTEM_OMP"])
env.AppendUnique(LIBS=['fftw3', 'fftw3_omp'])
conf = Configure(env)
if not conf.CheckLibWithHeader("Tamaas", "tamaas.hh", language="CXX"):
raise StopError(
"Tamaas not found in {}".format(os.path.abspath(tm_path)))
env['TAMAAS_PATH'] = os.path.abspath(tm_path)
env = conf.Finish()
def exists(env):
if 'TAMAAS_PATH' in env:
return True
return False
diff --git a/site_scons/version.py b/site_scons/version.py
index a5704db..d8e4a07 100644
--- a/site_scons/version.py
+++ b/site_scons/version.py
@@ -1,61 +1,62 @@
# -*- coding: utf-8 -*-
# @file
# @section LICENSE
#
# Copyright (©) 2016-19 EPFL (École Polytechnique Fédérale de Lausanne),
# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU Affero General Public License as published
# by the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Affero General Public License for more details.
#
# You should have received a copy of the GNU Affero General Public License
# along with this program. If not, see .
import subprocess
import base64
from functools import reduce
def git(*args):
return subprocess.check_output(["git"] + list(args))
def write_info_file(file_name, build_type="release"):
header = '#include "tamaas_info.hh"\n\n'
file_content = ''
branch = git("rev-parse", "--abbrev-ref", "HEAD")[:-1]
commit = git("rev-parse", branch)[:-1]
file_content += \
- 'std::string tamaas::TamaasInfo::build_type = "{}";\n'.format(build_type)
+ 'std::string tamaas::TamaasInfo::build_type = "{}";\n'.format(
+ build_type)
file_content += \
'std::string tamaas::TamaasInfo::branch = "{}";\n'.format(branch)
file_content += \
'std::string tamaas::TamaasInfo::commit = "{}";\n'.format(commit)
diff = git("diff")
diff += "|".encode('ascii') + git("diff", "--cached")
file_content += \
'std::string tamaas::TamaasInfo::diff = "{}";\n'.format(
base64.b64encode(diff))
remotes = git('remote', '-v')
remotes_strings = filter(lambda x: x != '', remotes.splitlines())
remotes_string = reduce(lambda x, y: x + '"{}\\n"\n'.format(y),
remotes_strings, "")
file_content += \
'std::string tamaas::TamaasInfo::remotes = {};\n'.format(remotes_string)
with open(file_name, 'w') as file:
file.write(header + file_content)
diff --git a/tests/fftfreq.py b/tests/fftfreq.py
index bfa61c3..71aaebf 100644
--- a/tests/fftfreq.py
+++ b/tests/fftfreq.py
@@ -1,54 +1,53 @@
# -*- coding: utf-8 -*-
# @file
# @section LICENSE
#
# Copyright (©) 2016-19 EPFL (École Polytechnique Fédérale de Lausanne),
# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU Affero General Public License as published
# by the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Affero General Public License for more details.
#
# You should have received a copy of the GNU Affero General Public License
# along with this program. If not, see .
from __future__ import division
__author__ = "Lucas Frérot"
import numpy as np
from numpy.fft import fftfreq
-from tamaas import dtype
def frequencies1D(frequencies):
n = len(frequencies)
frequencies[:] = n * fftfreq(n)
def assign2DFrequencies(qx, qy, frequencies):
qx, qy = np.meshgrid(qx, qy, indexing='ij')
frequencies[:, :, 0] = qx
frequencies[:, :, 1] = qy
def frequencies2D(frequencies):
n = frequencies.shape
qx = fftfreq(n[0]) * n[0]
qy = fftfreq(n[1]) * n[1]
assign2DFrequencies(qx, qy, frequencies)
def hfrequencies2D(frequencies):
n = frequencies.shape
qx = fftfreq(n[0]) * n[0]
qy = np.arange(n[1])
assign2DFrequencies(qx, qy, frequencies)
diff --git a/tests/test_autocorrelation.py b/tests/test_autocorrelation.py
index 6998446..a520be3 100644
--- a/tests/test_autocorrelation.py
+++ b/tests/test_autocorrelation.py
@@ -1,101 +1,101 @@
# -*- coding: utf-8 -*-
# @file
# @section LICENSE
#
# Copyright (©) 2016-19 EPFL (École Polytechnique Fédérale de Lausanne),
# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU Affero General Public License as published
# by the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Affero General Public License for more details.
#
# You should have received a copy of the GNU Affero General Public License
# along with this program. If not, see .
from __future__ import print_function
import tamaas as tm
import numpy as np
from scipy.special import jv
import scipy.integrate as integrate
def spectrum(r, kl, kr, ks, C, H):
if r > ks:
return 0
elif r > kr:
return C * (r/float(kr))**(-2.*(H+1))
elif r > kl:
return C
else:
return 0
# n is normalization parameter
def integrate_bessel(x, spectrum, n):
result = 2 * np.pi * \
integrate.quad(lambda y: y*jv(0, y*x) * spectrum(n*y), 0, np.inf)[0]
return result
def test_autocorrelation():
tm.initialize()
grid_size = 256
lambda_l = 8.
lambda_r = 8.
lambda_s = 1.
L = 4 * lambda_l
kl = int(L / lambda_l)
kr = int(L / lambda_r)
ks = int(L / lambda_s)
# Computing ACF from integral with Bessel function
- l = grid_size // 2 # points on curve
- x = np.linspace(0, L/2, l)
- bessel_acf = np.zeros(l)
+ n = grid_size // 2 # points on curve
+ x = np.linspace(0, L/2, n)
+ bessel_acf = np.zeros(n)
integrate_bessel_v = np.vectorize(integrate_bessel, otypes=[tm.dtype])
bessel_acf = (L / (2*np.pi))**2 *\
integrate_bessel_v(x, lambda x: spectrum(x,
kl, kr, ks,
1., 0.8), L / (2*np.pi))
nsamples = 1000
etages = np.zeros((grid_size, grid_size, nsamples))
# Computing ACF with Hu & Tonder algo reimplemented
for i in range(nsamples):
SG = tm.SurfaceGeneratorFilterFFT() # generate surface
SG.getGridSize().assign(grid_size)
SG.getHurst().assign(0.8)
SG.getQ0().assign(kl)
SG.getQ1().assign(kr)
SG.getQ2().assign(ks)
SG.getRandomSeed().assign(i)
SG.Init()
a = SG.buildSurface()
etages[:, :, i] = a[:, :]
- cov = np.zeros(l)
+ cov = np.zeros(n)
- for i in np.arange(l, grid_size, 1):
- cov[i-l] = np.sum((etages[l, i, :]-np.mean(etages[l, i, :])) *
- (etages[l, l, :]-np.mean(etages[l, l, :])))/nsamples
+ for i in np.arange(n, grid_size, 1):
+ cov[i-n] = np.sum((etages[n, i, :]-np.mean(etages[n, i, :])) *
+ (etages[n, n, :]-np.mean(etages[n, n, :])))/nsamples
error = np.sum((cov-bessel_acf)*(cov-bessel_acf)) \
/ np.sum(bessel_acf*bessel_acf)
print(error)
assert error < 1e-1
if __name__ == "__main__":
test_autocorrelation()
diff --git a/tests/test_bem_grid.py b/tests/test_bem_grid.py
index 04bee0f..95d8379 100644
--- a/tests/test_bem_grid.py
+++ b/tests/test_bem_grid.py
@@ -1,113 +1,112 @@
# -*- coding: utf-8 -*-
# @file
# @section LICENSE
#
# Copyright (©) 2016-19 EPFL (École Polytechnique Fédérale de Lausanne),
# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU Affero General Public License as published
# by the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Affero General Public License for more details.
#
# You should have received a copy of the GNU Affero General Public License
# along with this program. If not, see .
from __future__ import print_function
import tamaas as tm
import numpy as np
from numpy.linalg import norm
-from numpy.fft import fft2, ifft2
def make_profile(f, profile, k):
x = np.linspace(0, 1, profile.shape[0], endpoint=False, dtype=tm.dtype)
y = np.linspace(0, 1, profile.shape[1], endpoint=False, dtype=tm.dtype)
x, y = np.meshgrid(x, y)
profile[:, :] = f(2*k*np.pi*y)
def sin_profile(profile, k):
make_profile(np.sin, profile, k)
def cos_profile(profile, k):
make_profile(np.cos, profile, k)
def test_bem_grid():
"""Checking that all the pipes are clean"""
tm.initialize(1)
size = 128
k = 1
E = 0.91
nu = 0.3
mu = E / (2*(1+nu))
pressure = np.zeros((size, size), dtype=tm.dtype)
surface = np.zeros((size, size), dtype=tm.dtype)
westergaard_disp = np.zeros((size, size), dtype=tm.dtype)
cos_profile(pressure, k)
cos_profile(surface, k)
bem = tm.BemGridPolonski(surface)
bem.setElasticity(E, nu)
bem.computeInfluence()
bem.getTractions()[:, :, 2] = pressure[:, :]
bem.computeDisplacementsFromTractions()
displacements = bem.getDisplacements()
# Constructing vertical displacement solution
cos_profile(westergaard_disp, k)
westergaard_disp *= (1-nu**2)/(E*np.pi*k)
westergaard_norm = (1-nu**2)/(E*np.pi*k)*np.sqrt(0.5) * size
# Computing error
error = norm(displacements[:, :, 2] - westergaard_disp)/westergaard_norm
assert error < 1e-10, "Error in normal displacement"
# Constructing horizontal displacement solution
sin_profile(westergaard_disp, k)
westergaard_disp *= (2*nu-1)/(4*np.pi*k*mu)
westergaard_norm = np.abs(2*nu-1)/(4*np.pi*k*mu)*np.sqrt(0.5)*size
error = norm(displacements[:, :, 1] - westergaard_disp)/westergaard_norm
assert error < 1e-10, "Error in horizontal displacement"
# Looking at normal vectors
normals = np.zeros((size, size, 3), dtype=tm.dtype)
normals_norm = np.zeros((size, size, 3), dtype=tm.dtype)
sin_profile(normals[:, :, 1], k)
normals[:, :, 1] *= 2*np.pi*k
normals[:, :, 2] = 1
x = np.linspace(0, 1, size, endpoint=False, dtype=tm.dtype)
y = np.linspace(0, 1, size, endpoint=False, dtype=tm.dtype)
x, y = np.meshgrid(x, y)
normals_norm[:, :, 0] = np.sqrt(1+4*np.pi**2*k**2*np.sin(2*np.pi*k*y)**2)
normals_norm[:, :, 1] = np.sqrt(1+4*np.pi**2*k**2*np.sin(2*np.pi*k*y)**2)
normals_norm[:, :, 2] = np.sqrt(1+4*np.pi**2*k**2*np.sin(2*np.pi*k*y)**2)
normals /= normals_norm
bem.computeSurfaceNormals()
bem_normals = bem.getSurfaceNormals()
error = norm(bem_normals - normals)
assert error < 1e-9, "Error in surface normals"
tm.finalize()
if __name__ == "__main__":
test_bem_grid()
diff --git a/tests/test_flood_fill.py b/tests/test_flood_fill.py
index a414f63..8a685e7 100644
--- a/tests/test_flood_fill.py
+++ b/tests/test_flood_fill.py
@@ -1,71 +1,71 @@
# -*- coding: utf-8 -*-
# @file
# @section LICENSE
#
# Copyright (©) 2016-19 EPFL (École Polytechnique Fédérale de Lausanne),
# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU Affero General Public License as published
# by the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Affero General Public License for more details.
#
# You should have received a copy of the GNU Affero General Public License
# along with this program. If not, see .
from __future__ import print_function, division
import tamaas as tm
import numpy as np
def test_flood_fill2d(tamaas_fixture):
field = np.zeros((18, 18), np.bool)
# Cluster of permeter 18 area 13
field[2:4, 3:6] = True
- field[4 , 4:6] = True
- field[5 , 5 ] = True
- field[3:5, 6:8] = True
+ field[4 , 4:6] = True # noqa
+ field[5 , 5 ] = True # noqa
+ field[3:5, 6:8] = True # noqa
# Cluster of perimeter 8 area 4
field[14:16, 3:5] = True
# Cluster of perimeter 20 area 11
- field[9:11, 8:11 ] = True
- field[7:9 , 9 ] = True
- field[10 , 11:14] = True
+ field[9:11, 8:11 ] = True # noqa
+ field[7:9 , 9 ] = True # noqa
+ field[10 , 11:14] = True # noqa
# Cluster of perimeter 18 area 9
field[3:5, 14:16] = True
field[5:10, 15] = True
clusters = tm.FloodFill.getClusters(field, False)
# Dummy class for clusters
class dummy:
def __init__(self, area, perimeter):
self.area = area
self.perimeter = perimeter
# Solution
ref = [dummy(13, 18), dummy(9, 18), dummy(11, 20), dummy(4, 8)]
for x, y in zip(clusters, ref):
assert x.getArea() == y.area
assert x.getPerimeter() == y.perimeter
def test_flood_fill3d(tamaas_fixture):
field = np.zeros((5, 5, 5), np.bool)
field[2:4, 3, 1:3] = True
clusters = tm.FloodFill.getVolumes(field, False)
assert len(clusters) == 1
assert clusters[0].getArea() == 4
diff --git a/tests/test_gtest.py b/tests/test_gtest.py
index ab1d5ec..5fdcf15 100644
--- a/tests/test_gtest.py
+++ b/tests/test_gtest.py
@@ -1,34 +1,35 @@
# -*- coding: utf-8 -*-
# @file
# @section LICENSE
#
# Copyright (©) 2016-19 EPFL (École Polytechnique Fédérale de Lausanne),
# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU Affero General Public License as published
# by the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Affero General Public License for more details.
#
# You should have received a copy of the GNU Affero General Public License
# along with this program. If not, see .
import subprocess
import os
def test_gtest():
me_path = os.path.realpath(__file__)
current_dir = os.path.dirname(me_path)
my_env = os.environ.copy()
my_env['PYTHONPATH'] = current_dir + ":" + os.getenv('PYTHONPATH', "")
subprocess.check_call(
[os.path.join(current_dir, 'test_gtest_all'),
- '--gtest_output=xml:' + os.path.join(current_dir, 'gtest_results.xml')],
+ '--gtest_output=xml:' +
+ os.path.join(current_dir, 'gtest_results.xml')],
env=my_env
)
diff --git a/tests/test_hertz_disp.py b/tests/test_hertz_disp.py
index d984c6d..aac1d77 100644
--- a/tests/test_hertz_disp.py
+++ b/tests/test_hertz_disp.py
@@ -1,130 +1,136 @@
# -*- coding: utf-8 -*-
# @file
# @section LICENSE
#
# Copyright (©) 2016-19 EPFL (École Polytechnique Fédérale de Lausanne),
# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU Affero General Public License as published
# by the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Affero General Public License for more details.
#
# You should have received a copy of the GNU Affero General Public License
# along with this program. If not, see .
from __future__ import print_function, division
import tamaas as tm
import numpy as np
def constructHertzProfile(size, curvature):
radius = 1. / curvature
x = np.linspace(-0.5, 0.5, size, dtype=tm.dtype)
y = np.linspace(-0.5, 0.5, size, dtype=tm.dtype)
x, y = np.meshgrid(x, y)
surface = np.sqrt(radius**2 - x**2 - y**2)
surface -= surface.min()
return surface.copy()
def computeHertzDisplacement(e_star, contact_size, max_pressure, size):
x = np.linspace(-0.5, 0.5, size)
y = np.linspace(-0.5, 0.5, size)
x, y = np.meshgrid(x, y)
- disp = np.pi * max_pressure / (4 * contact_size * e_star) * (2 * contact_size**2 -
- (x**2 + y**2))
+ disp = np.pi * max_pressure / (4 * contact_size * e_star) \
+ * (2 * contact_size**2 - (x**2 + y**2))
return disp.copy()
def test_hertz_disp():
tm.initialize()
grid_size = 512
curvature = 0.5
effective_modulus = 1.
load = 0.12
surface = constructHertzProfile(grid_size, curvature)
model = tm.ModelFactory.createModel(tm.model_type_basic_2d,
[1., 1.],
[grid_size, grid_size])
model.setElasticity(1, 0)
solver = tm.PolonskyKeerRey(model, surface, 1e-12,
tm.PolonskyKeerRey.gap,
tm.PolonskyKeerRey.gap)
solver.solve(load - surface.mean())
tractions = model.getTraction()
true_displacements = model.getDisplacement()
true_gap = true_displacements - surface
pressure = np.mean(tractions)
contact_points = np.where(1.0e-10 >= true_gap)
contact_area = (np.size(contact_points)/2)/float(grid_size*grid_size)
print('The contact area computed with the gap map is '+str(contact_area))
- hertz_contact_size = (3 * pressure / (4 * curvature * effective_modulus))**(1. / 3.)
+ hertz_contact_size = (3 * pressure
+ / (4 * curvature * effective_modulus))**(1. / 3.)
print('Exact Hertz contact radius is '+str(hertz_contact_size))
hertz_area = np.pi * hertz_contact_size**2
print('Exact Hertz contact area is '+str(hertz_area))
area_error = np.abs(hertz_area - contact_area) / hertz_area
print("Area error: {}".format(area_error))
# Testing maximum pressure
max_pressure = tractions.max()
- hertz_max_pressure = (6 * pressure * effective_modulus**2 * curvature ** 2)**(1. / 3.) / np.pi
- pressure_error = np.abs(hertz_max_pressure - max_pressure) / hertz_max_pressure
+ hertz_max_pressure = (6 * pressure * effective_modulus**2
+ * curvature ** 2)**(1. / 3.) / np.pi
+ pressure_error = np.abs(hertz_max_pressure - max_pressure) \
+ / hertz_max_pressure
print('Exact Hertz max pressure is '+str(hertz_max_pressure))
print('max pressure is '+str(np.max(tractions)))
print("Max pressure error: {}".format(pressure_error))
# Testing displacements
hertz_displacements = computeHertzDisplacement(effective_modulus,
hertz_contact_size,
hertz_max_pressure,
grid_size)
# Selecing only the points that are in contact
contact_indexes = [(i, j)
for i in range(grid_size) for j in range(grid_size)
if true_gap[i, j] < 1e-12]
# Displacements of bem are centered around the mean of the whole surface
# and Hertz displacements are not centered, so we need to compute mean
# on the contact zone for both arrays
bem_mean = 0.
hertz_mean = 0.
for index in contact_indexes:
bem_mean += true_displacements[index]
hertz_mean += hertz_displacements[index]
bem_mean /= len(contact_indexes)
hertz_mean /= len(contact_indexes)
# Correction applied when computing error
correction = hertz_mean - bem_mean
# Computation of error of displacement in contact zone
error = 0.
hertz_norm = 0.
for index in contact_indexes:
- error += (hertz_displacements[index] - true_displacements[index] - correction)**2
+ error += (hertz_displacements[index] -
+ true_displacements[index] - correction)**2
hertz_norm += (hertz_displacements[index] - hertz_mean)**2
displacement_error = np.sqrt(error / hertz_norm)
print("Displacement error (in contact zone): {}".format(displacement_error))
- assert area_error < 1e-2 and pressure_error < 1e-2 and displacement_error < 2e-3
+ assert area_error < 1e-2 \
+ and pressure_error < 1e-2 \
+ and displacement_error < 2e-3
tm.finalize()
if __name__ == "__main__":
test_hertz_disp()
diff --git a/tests/test_hertz_kato.py b/tests/test_hertz_kato.py
index b224ad6..cb70274 100644
--- a/tests/test_hertz_kato.py
+++ b/tests/test_hertz_kato.py
@@ -1,117 +1,123 @@
# -*- coding: utf-8 -*-
# @file
# @section LICENSE
#
# Copyright (©) 2016-19 EPFL (École Polytechnique Fédérale de Lausanne),
# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU Affero General Public License as published
# by the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Affero General Public License for more details.
#
# You should have received a copy of the GNU Affero General Public License
# along with this program. If not, see .
from __future__ import division, print_function
import numpy as np
import tamaas as tm
def constructHertzProfile(size, curvature):
radius = 1. / curvature
x = np.linspace(-0.5, 0.5, size, dtype=tm.dtype)
y = np.linspace(-0.5, 0.5, size, dtype=tm.dtype)
x, y = np.meshgrid(x, y)
surface = np.sqrt(radius**2 - x**2 - y**2)
surface -= surface.mean()
return surface.copy()
def computeHertzDisplacement(e_star, contact_size, max_pressure, size):
x = np.linspace(-0.5, 0.5, size, dtype=tm.dtype)
y = np.linspace(-0.5, 0.5, size, dtype=tm.dtype)
x, y = np.meshgrid(x, y)
- disp = np.pi * max_pressure / (4 * contact_size * e_star) * (2 * contact_size**2 -
- (x**2 + y**2))
+ disp = np.pi * max_pressure / (4 * contact_size * e_star) \
+ * (2 * contact_size**2 - (x**2 + y**2))
return disp.copy()
def test_hertz_kato():
tm.initialize()
grid_size = 512
curvature = 0.1
effective_modulus = 1.
load = 0.0001
surface = constructHertzProfile(grid_size, curvature)
bem = tm.BemKato(surface)
bem.setEffectiveModulus(effective_modulus)
bem.setMaxIterations(100000)
bem.computeEquilibrium(1e-7, load)
tractions = bem.getTractions()
displacements = bem.getDisplacements()
# Testing contact area against Hertz solution for solids of revolution
contact_area = tm.SurfaceStatistics.computeContactArea(tractions)
- hertz_contact_size = (3 * load / (4 * curvature * effective_modulus))**(1. / 3.)
+ hertz_contact_size = (3 * load / (4 * curvature
+ * effective_modulus))**(1. / 3.)
hertz_area = np.pi * hertz_contact_size**2
area_error = np.abs(hertz_area - contact_area) / hertz_area
print("Area error: {}".format(area_error))
# Testing maximum pressure
max_pressure = tractions.max()
- hertz_max_pressure = (6 * load * effective_modulus**2 * curvature ** 2)**(1. / 3.) / np.pi
- pressure_error = np.abs(hertz_max_pressure - max_pressure) / hertz_max_pressure
+ hertz_max_pressure = (6 * load * effective_modulus**2
+ * curvature ** 2)**(1. / 3.) / np.pi
+ pressure_error = np.abs(hertz_max_pressure - max_pressure) \
+ / hertz_max_pressure
print("Max pressure error: {}".format(pressure_error))
# Testing displacements
hertz_displacements = computeHertzDisplacement(effective_modulus,
hertz_contact_size,
hertz_max_pressure,
grid_size)
# Selecing only the points that are in contact
contact_indexes = [(i, j)
for i in range(grid_size) for j in range(grid_size)
if tractions[i, j] > 0]
# Displacements of bem are centered around the mean of the whole surface
# and Hertz displacements are not centered, so we need to compute mean
# on the contact zone for both arrays
bem_mean = 0.
hertz_mean = 0.
for index in contact_indexes:
bem_mean += displacements[index]
hertz_mean += hertz_displacements[index]
bem_mean /= len(contact_indexes)
hertz_mean /= len(contact_indexes)
# Correction applied when computing error
correction = hertz_mean - bem_mean
# Computation of error of displacement in contact zone
error = 0.
hertz_norm = 0.
for index in contact_indexes:
- error += (hertz_displacements[index] - displacements[index] - correction)**2
+ error += (hertz_displacements[index]
+ - displacements[index] - correction)**2
hertz_norm += (hertz_displacements[index] - hertz_mean)**2
displacement_error = np.sqrt(error / hertz_norm)
print("Displacement error (in contact zone): {}".format(displacement_error))
- assert area_error < 6e-3 and pressure_error < 4e-3 and displacement_error < 3e-4
+ assert area_error < 6e-3 \
+ and pressure_error < 4e-3 \
+ and displacement_error < 3e-4
tm.finalize()
if __name__ == "__main__":
test_hertz_kato()
diff --git a/tests/test_tangential.py b/tests/test_tangential.py
index 791bdfa..760b3c2 100644
--- a/tests/test_tangential.py
+++ b/tests/test_tangential.py
@@ -1,102 +1,106 @@
# -*- coding: utf-8 -*-
# @file
# @section LICENSE
#
# Copyright (©) 2016-19 EPFL (École Polytechnique Fédérale de Lausanne),
# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU Affero General Public License as published
# by the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Affero General Public License for more details.
#
# You should have received a copy of the GNU Affero General Public License
# along with this program. If not, see .
from __future__ import division, print_function
import numpy as np
from numpy.linalg import norm
import tamaas as tm
def test_pressure():
tm.initialize()
E = 1.0
nu = 0.5
mu = 0.5
n = 81
L = 1.0
p_target = np.array([0.5e-4, 0, 2e-4], dtype=tm.dtype)
- g_target = np.array([-0.000223, 0, 9.99911], dtype=tm.dtype)
+ # g_target = np.array([-0.000223, 0, 9.99911], dtype=tm.dtype)
x = np.linspace(-L/2, L/2, n, dtype=tm.dtype)
y = np.copy(x)
x, y = np.meshgrid(x, y, indexing="ij")
r_sqrd = (x**2 + y**2)
# Spherical contact
R = 10.0
surface = R * np.sqrt((r_sqrd < R**2) * (1 - r_sqrd / R**2))
# Creating model
- model = tm.ModelFactory.createModel(tm.model_type_surface_2d, [L, L], [n, n])
+ model = tm.ModelFactory.createModel(tm.model_type_surface_2d,
+ [L, L], [n, n])
model.setElasticity(E, nu)
p = model.getTraction()
# Theoretical solution
Estar = E / (1.0 - nu**2)
Fn = p_target[2] * L**2
Fx = p_target[0] * L**2
d_theory = np.power(3 / 4 * Fn / Estar / np.sqrt(R), 2/3)
a_theory = np.sqrt(R * d_theory)
c_theory = a_theory * (1 - Fx / (mu * Fn)) ** (1/3)
p0_theory = np.power(6 * Fn * Estar**2 / np.pi**3 / R**2, 1/3)
t1_theory = mu * p0_theory
t2_theory = t1_theory * c_theory / a_theory
- # p_theory = p0_theory * np.sqrt((r_sqrd < a_theory**2) * (1 - r_sqrd / a_theory**2))
- t_theory = t1_theory * np.sqrt((r_sqrd < a_theory**2) * (1 - r_sqrd / a_theory**2))
- t_theory = t_theory - t2_theory * np.sqrt((r_sqrd < c_theory**2) * (1 - r_sqrd / c_theory**2))
+ # p_theory = p0_theory * np.sqrt((r_sqrd < a_theory**2)
+ # * (1 - r_sqrd / a_theory**2))
+ t_theory = t1_theory * np.sqrt((r_sqrd < a_theory**2)
+ * (1 - r_sqrd / a_theory**2))
+ t_theory = t_theory - t2_theory * np.sqrt((r_sqrd < c_theory**2)
+ * (1 - r_sqrd / c_theory**2))
def assert_error(err):
error = norm(p[..., 0] - t_theory) / norm(t_theory)
print(error)
assert error < err
# Test Kato solver
solver = tm.Kato(model, surface, 1e-12, mu)
solver.setMaxIterations(1000)
solver.solve(p_target, 100)
assert_error(4e-2)
# solver.solveRelaxed(g_target)
# assert_error(1e-2)
# # Test BeckTeboulle solver
# solver = tm.BeckTeboulle(model, surface, 1e-12, mu)
# solver.solve(g_target)
# assert_error(1e-2)
# Test Condat solver
solver = tm.Condat(model, surface, 1e-12, mu)
- solver.setMaxIterations(5000) # or 10000
+ solver.setMaxIterations(5000) # or 10000
solver.solve(p_target)
- assert_error(7e-2) # 4e-2 for 10000 iterations
+ assert_error(7e-2) # 4e-2 for 10000 iterations
# Test tangential Polonsky Keer solver
solver = tm.PolonskyKeerTan(model, surface, 1e-12, mu)
solver.setMaxIterations(1000)
solver.solve(p_target)
assert_error(4e-2)
tm.finalize()
if __name__ == "__main__":
test_pressure()
diff --git a/tests/test_voigt.py b/tests/test_voigt.py
index e1cdfd5..8abfb5d 100644
--- a/tests/test_voigt.py
+++ b/tests/test_voigt.py
@@ -1,34 +1,34 @@
# -*- coding: utf-8 -*-
# @file
# @section LICENSE
#
# Copyright (©) 2016-19 EPFL (École Polytechnique Fédérale de Lausanne),
# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU Affero General Public License as published
# by the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Affero General Public License for more details.
#
# You should have received a copy of the GNU Affero General Public License
# along with this program. If not, see .
import tamaas as tm
import numpy as np
def test_voigt(tamaas_fixture):
shape = [3, 3, 3]
stress = np.ones(shape + [9], dtype=tm.dtype)
- #voigt = np.zeros(shape + [6])
+ # voigt = np.zeros(shape + [6])
voigt = tm.to_voigt(stress)
sol = np.ones(shape + [6], dtype=tm.dtype)
sol[:, :, :, 3:] = np.sqrt(2)
print(voigt)
print(sol)
assert np.all(np.abs(sol - voigt) < 1e-15)
diff --git a/tests/test_westergaard.py b/tests/test_westergaard.py
index 900ad4c..e66ba23 100644
--- a/tests/test_westergaard.py
+++ b/tests/test_westergaard.py
@@ -1,80 +1,78 @@
# -*- coding: utf-8 -*-
# @file
# @section LICENSE
#
# Copyright (©) 2016-19 EPFL (École Polytechnique Fédérale de Lausanne),
# Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides)
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU Affero General Public License as published
# by the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Affero General Public License for more details.
#
# You should have received a copy of the GNU Affero General Public License
# along with this program. If not, see .
from __future__ import division, print_function
-import sys
-import numpy as np
from numpy.linalg import norm
import tamaas as tm
import pytest
@pytest.fixture(scope="module",
params=[tm.PolonskyKeerRey.pressure,
tm.PolonskyKeerRey.gap])
def pkr(westergaard, request):
loads = {
tm.PolonskyKeerRey.pressure: westergaard.load,
tm.PolonskyKeerRey.gap: 0.06697415181446396
}
model = tm.ModelFactory.createModel(tm.model_type.basic_1d, [1.],
[westergaard.n])
model.setElasticity(westergaard.e_star, 0)
solver = tm.PolonskyKeerRey(model, westergaard.surface, 1e-12,
request.param,
request.param)
solver.setMaxIterations(5000)
solver.solve(loads[request.param])
return model, westergaard
def test_energy(pkr):
model, sol = pkr
traction, displacement = model.getTraction(), model.getDisplacement()
error = norm((traction - sol.pressure)*(displacement - sol.displacement))
error /= norm(sol.pressure * sol.displacement)
assert error < 4e-6
def test_pkr_multi(westergaard):
model_basic = tm.ModelFactory.createModel(tm.model_type.basic_1d, [1.],
[westergaard.n])
model_volume = tm.ModelFactory.createModel(tm.model_type.volume_1d,
[1., 1.],
[20, westergaard.n])
pressures = {}
for model in [model_basic, model_volume]:
print(model)
solver = tm.PolonskyKeerRey(model, westergaard.surface, 1e-12,
tm.PolonskyKeerRey.pressure,
tm.PolonskyKeerRey.pressure)
solver.solve(westergaard.load)
pressures[model] = model.getTraction()
error = norm(pressures[model_basic] - pressures[model_volume][:, 1])\
/ westergaard.load
assert error < 1e-16
if __name__ == "__main__":
print("Test is meant to run with pytest")