From 559d9c09dbcfa9952160d208dd5f4e271228a061 Mon Sep 17 00:00:00 2001 From: Mark Wieczorek Date: Sun, 31 Mar 2024 00:03:05 +0100 Subject: [PATCH 1/7] Add area and area_equivalent_radius for sphere and ellipsoid --- boule/_ellipsoid.py | 34 ++++++++++++++++++++++++++++++++++ boule/_sphere.py | 22 ++++++++++++++++++++++ 2 files changed, 56 insertions(+) diff --git a/boule/_ellipsoid.py b/boule/_ellipsoid.py index c2f9a820..77ea1853 100644 --- a/boule/_ellipsoid.py +++ b/boule/_ellipsoid.py @@ -103,6 +103,10 @@ class Ellipsoid: 6371008.7714 m >>> print(f"{ellipsoid.volume * 1e-9:.5e} km³") 1.08321e+12 km³ + >>> print(f"{ellipsoid.area:.10e} m²") + 5.1006562172e+14 m² + >>> print(f"{ellipsoid.area_equivalent_radius:0.4f} m") + 6371007.1809 m >>> print(f"{ellipsoid.gravity_equator:.10f} m/s²") 9.7803253359 m/s² >>> print(f"{ellipsoid.gravity_pole:.10f} m/s²") @@ -219,6 +223,27 @@ def mean_radius(self): """ return 1 / 3 * (2 * self.semimajor_axis + self.semiminor_axis) + @property + def area(self): + r""" + The area of the ellipsoid. + Definition: :math:`A = 2 \pi a^2 \left(1 + \dfrac{b^2}{e a^2} + \arctanh{e} \right)`. + Units: :math:`m^2`. + """ + # see https://en.wikipedia.org/wiki/Ellipsoid#Surface_area + return ( + 2 + * np.pi + * self.semimajor_axis**2 + * ( + 1 + + (self.semiminor_axis / self.semimajor_axis) ** 2 + / self.first_eccentricity + * np.arctanh(self.first_eccentricity) + ) + ) + @property def volume(self): r""" @@ -228,6 +253,15 @@ def volume(self): """ return (4 / 3 * np.pi) * self.semimajor_axis**2 * self.semiminor_axis + @property + def area_equivalent_radius(self): + r""" + The area equivalent radius of the ellipsoid. + Definition: :math:`R_2 = \sqrt{A / (4 \pi)}`. + Units: :math:`m`. + """ + return np.sqrt(self.area / 4 / np.pi) + @property def _emm(self): "Auxiliary quantity used to calculate gravity at the pole and equator" diff --git a/boule/_sphere.py b/boule/_sphere.py index 462a1b04..9398460a 100644 --- a/boule/_sphere.py +++ b/boule/_sphere.py @@ -99,6 +99,10 @@ class Sphere: 0 >>> print(f"{sphere.volume:.10f} m³") 4.1887902048 m³ + >>> print(f"{sphere.area:.10f} m²") + 12.5663706144 m² + >>> print(sphere.area_equivalent_radius) + 1 """ @@ -176,6 +180,15 @@ def first_eccentricity(self): """ return 0 + @property + def area(self): + r""" + The area of the sphere. + Definition: :math:`A = 4 \pi r^2`. + Units: :math:`m^2`. + """ + return 4 * np.pi * self.radius**2 + @property def volume(self): r""" @@ -185,6 +198,15 @@ def volume(self): """ return (4 / 3 * np.pi) * self.radius**3 + @property + def area_equivalent_radius(self): + """ + The area equivalent radius of the sphere is equal to its radius. + Definition: :math:`R_2 = R`. + Units: :math:`m`. + """ + return self.radius + def normal_gravity(self, latitude, height, si_units=False): r""" Normal gravity of the sphere at the given latitude and height. From 42918247350e323e87c8569c965c773fb5d2574d Mon Sep 17 00:00:00 2001 From: Mark Wieczorek Date: Sun, 31 Mar 2024 19:31:59 +0200 Subject: [PATCH 2/7] Add area and R2 to TriaxialEllipsoid --- boule/_triaxialellipsoid.py | 33 +++++++++++++++++++++++++++++++++ environment.yml | 1 + 2 files changed, 34 insertions(+) diff --git a/boule/_triaxialellipsoid.py b/boule/_triaxialellipsoid.py index d3fe2d6e..2b91e0df 100644 --- a/boule/_triaxialellipsoid.py +++ b/boule/_triaxialellipsoid.py @@ -11,6 +11,7 @@ import attr import numpy as np +from scipy.special import elliprd, elliprf # Don't let ellipsoid parameters be changed to avoid messing up calculations @@ -98,6 +99,10 @@ class TriaxialEllipsoid: >>> print(f"{ellipsoid.mean_radius:.0f} m") 262700 m + >>> print(f"{ellipsoid.area:.10e} m²") + 8.6562393883e+11 m² + >>> print(f"{ellipsoid.area_equivalent_radius:0.0f} m") + 262458 m >>> print(f"{ellipsoid.volume * 1e-9:.0f} km³") 74573626 km³ @@ -171,6 +176,25 @@ def mean_radius(self): """ return (self.semimajor_axis + self.semimedium_axis + self.semiminor_axis) / 3 + @property + def area(self): + r""" + The area of the ellipsoid. + Definition: :math:`A`. + Units: :math:`m^2`. + """ + # see https://en.wikipedia.org/wiki/Ellipsoid#Surface_area + a = self.semimajor_axis + b = self.semimedium_axis + c = self.semiminor_axis + return 2 * np.pi * c**2 + 2 * np.pi * a * b * ( + elliprf((c / a) ** 2, (c / b) ** 2, 1) + - (1 / 3) + * (1 - (c / a) ** 2) + * (1 - (c / b) ** 2) + * elliprd((c / a) ** 2, (c / b) ** 2, 1) + ) + @property def volume(self): r""" @@ -185,6 +209,15 @@ def volume(self): * self.semiminor_axis ) + @property + def area_equivalent_radius(self): + r""" + The area equivalent radius of the ellipsoid. + Definition: :math:`R_2 = \sqrt{A / (4 \pi)}`. + Units: :math:`m`. + """ + return np.sqrt(self.area / 4 / np.pi) + @property def equatorial_flattening(self): r""" diff --git a/environment.yml b/environment.yml index e821a715..b73a90dd 100644 --- a/environment.yml +++ b/environment.yml @@ -9,6 +9,7 @@ dependencies: # Run - numpy - attrs + - scipy # Build - build - twine From b249a33ac8785243b57d60c68f30dce12ea5e884 Mon Sep 17 00:00:00 2001 From: Mark Wieczorek Date: Sun, 31 Mar 2024 19:48:22 +0200 Subject: [PATCH 3/7] Add scipy to requirements-test.txt --- env/requirements-test.txt | 1 + 1 file changed, 1 insertion(+) diff --git a/env/requirements-test.txt b/env/requirements-test.txt index 180030ae..beebd9c9 100644 --- a/env/requirements-test.txt +++ b/env/requirements-test.txt @@ -3,3 +3,4 @@ pytest pytest-cov coverage pymap3d>=2.9.0 +scipy From e6b03b89fe338a1951f1b2b3b1e80fbbe70f55d5 Mon Sep 17 00:00:00 2001 From: Mark Wieczorek Date: Fri, 12 Apr 2024 12:03:48 +0200 Subject: [PATCH 4/7] Update python tests to use 3.8 and 3.12 --- .github/workflows/test.yml | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/.github/workflows/test.yml b/.github/workflows/test.yml index 17efe862..e0e5d47a 100644 --- a/.github/workflows/test.yml +++ b/.github/workflows/test.yml @@ -38,11 +38,11 @@ jobs: fail-fast: false matrix: os: [ubuntu, macos, windows] - python: ["3.7", "3.11"] + python: ["3.8", "3.12"] include: - - python: "3.7" + - python: "3.8" dependencies: oldest - - python: "3.11" + - python: "3.12" dependencies: latest env: REQUIREMENTS: env/requirements-build.txt env/requirements-test.txt From b6c7af3bf5dd43ef52c1cde71f7c6b4c9dce153c Mon Sep 17 00:00:00 2001 From: Mark Wieczorek Date: Tue, 16 Apr 2024 17:28:57 +0200 Subject: [PATCH 5/7] Add name to AUTHORS.md --- AUTHORS.md | 1 + 1 file changed, 1 insertion(+) diff --git a/AUTHORS.md b/AUTHORS.md index 9c1bd9c1..fd96c705 100644 --- a/AUTHORS.md +++ b/AUTHORS.md @@ -10,3 +10,4 @@ order by last name) and are considered "The Boule Developers": * [Agustina Pesce](https://github.com/aguspesce) - Unaffiliated (ORCID: [0000-0002-5538-8845](https://orcid.org//0000-0002-5538-8845)) * [Santiago Soler](https://github.com/santisoler) - CONICET, Argentina; Instituto Geofísico Sismológico Volponi, Universidad Nacional de San Juan, Argentina (ORCID: [0000-0001-9202-5317](https://www.orcid.org/0000-0001-9202-5317)) * [Leonardo Uieda](https://github.com/leouieda) - Universidade de São Paulo, Brazil (ORCID: 0000-0001-6123-9515) +* [Mark Wieczorek](https://github.com/MarkWieczorek) - Institut de physique du globe de Paris, France (ORCID: [0000-0001-7007-4222](https://orcid.org/0000-0001-7007-4222)) From 968bbf28c3ea60e48c7f304b3dd7be518c8d57f4 Mon Sep 17 00:00:00 2001 From: Mark Wieczorek Date: Tue, 16 Apr 2024 22:46:53 +0200 Subject: [PATCH 6/7] Add minimum scipy version of 1.8.0 --- pyproject.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pyproject.toml b/pyproject.toml index 1a02857c..4e94f6f9 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -31,7 +31,7 @@ requires-python = ">=3.8" dependencies = [ "numpy>=1.19", "attrs>=19.3.0", - "scipy>=1.0.8", + "scipy>=1.8.0", ] [project.urls] From c8031291bc6383a3fc4560f47f4859b9896f4165 Mon Sep 17 00:00:00 2001 From: Leonardo Uieda Date: Wed, 17 Apr 2024 09:27:40 -0300 Subject: [PATCH 7/7] Use the R_G version of the area equation --- boule/_ellipsoid.py | 2 +- boule/_triaxialellipsoid.py | 23 ++++++++++++----------- 2 files changed, 13 insertions(+), 12 deletions(-) diff --git a/boule/_ellipsoid.py b/boule/_ellipsoid.py index 7842d478..922866f9 100644 --- a/boule/_ellipsoid.py +++ b/boule/_ellipsoid.py @@ -264,7 +264,7 @@ def area(self): r""" The area of the ellipsoid. Definition: :math:`A = 2 \pi a^2 \left(1 + \dfrac{b^2}{e a^2} - \arctanh{e} \right)`. + \text{arctanh}\,e \right)`. Units: :math:`m^2`. """ # see https://en.wikipedia.org/wiki/Ellipsoid#Surface_area diff --git a/boule/_triaxialellipsoid.py b/boule/_triaxialellipsoid.py index d1757ca8..58f712e1 100644 --- a/boule/_triaxialellipsoid.py +++ b/boule/_triaxialellipsoid.py @@ -11,7 +11,7 @@ import attr import numpy as np -from scipy.special import elliprd, elliprf +import scipy.special from ._constants import G @@ -231,19 +231,20 @@ def semiaxes_mean_radius(self): def area(self): r""" The area of the ellipsoid. - Definition: :math:`A`. + Definition: :math:`A = 3 V R_G(a^{-2}, b^{-2}, c^{-2})`, in which + :math:`R_G` is the completely-symmetric elliptic integral of the second + kind. Units: :math:`m^2`. """ # see https://en.wikipedia.org/wiki/Ellipsoid#Surface_area - a = self.semimajor_axis - b = self.semimedium_axis - c = self.semiminor_axis - return 2 * np.pi * c**2 + 2 * np.pi * a * b * ( - elliprf((c / a) ** 2, (c / b) ** 2, 1) - - (1 / 3) - * (1 - (c / a) ** 2) - * (1 - (c / b) ** 2) - * elliprd((c / a) ** 2, (c / b) ** 2, 1) + return ( + 3 + * self.volume + * scipy.special.elliprg( + 1 / self.semimajor_axis**2, + 1 / self.semimedium_axis**2, + 1 / self.semiminor_axis**2, + ) ) @property