Performance mode

.gitignore

@@ -4,6 +4,7 @@
 *.code-workspace
 .DS_Store
 docs/build/
+LocalPreferences.toml
 Manifest.toml
 tmp/
 .vscode/

NEWS.md

@@ -13,6 +13,10 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 
 - Define `GeneralPolytope` that represents general polytopes in 2 and 3 dimensions. Since PR[#1006](https://github.com/gridap/Gridap.jl/pull/1006).
 
+### Changed
+
+- Added local preferences for Gridap through `Preferences.jl`. The macro `@check` can now be activated/deactivated by using the local preference `execution_mode`. Since PR[#1014](https://github.com/gridap/Gridap.jl/pull/1014).
+
 ## [0.18.2] - 2024-05-02
 
 ### Fixed

Project.toml

@@ -21,6 +21,7 @@ LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 NLsolve = "2774e3e8-f4cf-5e23-947b-6d7e65073b56"
 NearestNeighbors = "b8a86587-4115-5ab1-83bc-aa920d37bbce"
 PolynomialBases = "c74db56a-226d-5e98-8bb0-a6049094aeea"
+Preferences = "21216c6a-2e73-6563-6e65-726566657250"
 QuadGK = "1fd47b50-473d-5c70-9696-f719f8f3bcdc"
 Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
 SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
@@ -47,6 +48,7 @@ LineSearches = "7.0.1"
 NLsolve = "4.3.0"
 NearestNeighbors = "0.4.8"
 PolynomialBases = "0.4.12"
+Preferences = "1.4"
 QuadGK = "2.3.1, 2.4"
 SparseMatricesCSR = "0.6.4"
 StaticArrays = "0.12.1, 1.0"

README.md

@@ -60,6 +60,7 @@ pkg> add Gridap
 ```
 
 ## Plugins
+
 - [GridapDistributed](https://github.com/gridap/GridapDistributed.jl) Distributed-memory extension of Gridap.
 - [GridapEmbedded](https://github.com/gridap/GridapEmbedded.jl) Embedded finite elements in Julia.
 - [GridapGmsh](https://github.com/gridap/GridapGmsh.jl) Generate a FE mesh with [GMSH](www.gmsh.info) and use it in Gridap.

src/Adaptivity/EdgeBasedRefinement.jl

@@ -1020,7 +1020,7 @@
 end
 
 function get_relabeled_connectivity(::PowellSabinRefinement,rr::RefinementRule,faces_gids)
-  @assert is_simplex(get_polytope(rr))
+  @check is_simplex(get_polytope(rr))
   conn = rr.ref_grid.grid.cell_node_ids
   gids = vcat(faces_gids...)
   new_data = lazy_map(Reindex(gids),conn.data)

src/Adaptivity/OldToNewFields.jl

@@ -18,7 +18,7 @@
 function OldToNewField(old_fields::AbstractArray{<:Fields.Field},
                        rrule::RefinementRule,
                        child_ids::AbstractVector{<:Integer})
-  @assert length(old_fields)==length(child_ids)                   
+  @check length(old_fields)==length(child_ids)                   
   if length(old_fields) == 1
     cell_map  = get_cell_map(rrule)[child_ids[1]]
     old_field = old_fields[1]
@@ -29,7 +29,7 @@
     refined_or_untouched_field = old_field∘cell_map
     return OldToNewField(RefinedOrUntouchedNewFieldType(),fine_to_coarse_field,refined_or_untouched_field)
   else 
-    @assert length(old_fields) <= num_subcells(rrule)
+    @check length(old_fields) <= num_subcells(rrule)
     fine_to_coarse_field = FineToCoarseField(old_fields,rrule,child_ids)
     cell_map = get_cell_map(rrule)[1]
     refined_or_untouched_field = old_fields[1]∘cell_map 

src/Arrays/Tables.jl

@@ -276,7 +276,7 @@ function append_tables_locally(offsets::NTuple, tables::NTuple)
 
   first_table, = tables
 
-  @assert all( map(length,tables) .== length(first_table) ) "All tables must have the same length"
+  @check all( map(length,tables) .== length(first_table) ) "All tables must have the same length"
   ndata = sum( (length(table.data) for table in tables) )
 
   T = get_data_eltype(first_table)

src/CellData/CellDofs.jl

@@ -35,7 +35,7 @@
   if strian === ttrian
     return a
   end
-  @assert is_change_possible(strian,ttrian) msg
+  @check is_change_possible(strian,ttrian) msg
   D = num_cell_dims(strian)
   sglue = get_glue(strian,Val(D))
   tglue = get_glue(ttrian,Val(D))
@@ -51,7 +51,7 @@
   if strian === ttrian
     return a
   end
-  @assert is_change_possible(strian,ttrian) msg
+  @check is_change_possible(strian,ttrian) msg
   D = num_cell_dims(strian)
   sglue = get_glue(strian,Val(D))
   tglue = get_glue(ttrian,Val(D))

src/CellData/CellFields.jl

@@ -156,7 +156,7 @@ function change_domain(a::CellField,strian::Triangulation,::ReferenceDomain,ttri
   if strian === ttrian
     return a
   end
-  @assert is_change_possible(strian,ttrian) msg
+  @check is_change_possible(strian,ttrian) msg
   D = num_cell_dims(strian)
   sglue = get_glue(strian,Val(D))
   tglue = get_glue(ttrian,Val(D))
@@ -172,7 +172,7 @@ function change_domain(a::CellField,strian::Triangulation,::PhysicalDomain,ttria
   if strian === ttrian
     return a
   end
-  @assert is_change_possible(strian,ttrian) msg
+  @check is_change_possible(strian,ttrian) msg
   D = num_cell_dims(strian)
   sglue = get_glue(strian,Val(D))
   tglue = get_glue(ttrian,Val(D))

src/FESpaces/CLagrangianFESpaces.jl

@@ -194,7 +194,7 @@ function _generate_node_to_dof_glue_component_major(
   nfree_dofs = 0
   ndiri_dofs = 0
   ncomps = length(z)
-  @assert length(testitem(tag_to_masks)) == ncomps
+  @check length(testitem(tag_to_masks)) == ncomps
   for (node,tag) in enumerate(node_to_tag)
     if tag == UNSET
       nfree_dofs += ncomps

src/FESpaces/FESpaceInterface.jl

@@ -23,7 +23,7 @@ function get_cell_fe_data(fun,f,ttrian)
   if strian === ttrian
     return sface_to_data
   end
-  @assert is_change_possible(strian,ttrian)
+  @check is_change_possible(strian,ttrian)
   D = num_cell_dims(strian)
   sglue = get_glue(strian,Val(D))
   tglue = get_glue(ttrian,Val(D))

src/Helpers/Helpers.jl

@@ -9,6 +9,10 @@ $(EXPORTS)
 module Helpers
 using DocStringExtensions
 
+@static if VERSION >= v"1.6"
+  using Preferences
+end
+
 import Base: +, -, *, /, transpose, adjoint
 import LinearAlgebra: cross, tr, dot
 
@@ -21,8 +25,11 @@ export tfill
 export get_val_parameter
 export first_and_tail
 export GridapType
+export set_debug_mode, set_performance_mode
 #export operate
 
+include("Preferences.jl")
+
 include("Macros.jl")
 
 include("HelperFunctions.jl")

src/Helpers/Macros.jl

@@ -53,11 +53,16 @@ end
   @check condition "Error message"
 
 Macro used to make sure that condition is fulfilled, like `@assert`
-but the check gets deactivated when running Julia with --boundscheck=no
+but the check gets deactivated when running Gridap in performance mode.
 """
 macro check(test,msg="A check failed")
-  quote
-    @boundscheck @assert $(esc(test)) $(esc(msg))
+  @static if execution_mode == "debug"
+    quote
+      @assert $(esc(test)) $(esc(msg))
+    end
+  else
+    quote
+      nothing
+    end
   end
 end
-

src/Helpers/Preferences.jl

@@ -0,0 +1,45 @@
+
+"""
+    set_execution_mode(new_mode::String)
+
+Sets the execution mode to either "debug" or "performance", which controls 
+the behavior of the @check macro within the Gridap package.
+
+- Debug mode (default): The @check macro will be active, which activates consistency checks 
+  within the library. This mode is recommended for development and debugging purposes.
+
+- Performance mode: The @check macro will be deactivated. This mode is recommended for 
+  production runs, where no errors are expected.
+
+Pre-defined functions [`set_debug_mode`](@ref) and [`set_performance_mode`](@ref) are also available.
+Feature only available in Julia 1.6 and later due to restrictions from `Preferences.jl`.
+"""
+function set_execution_mode(new_mode::String)
+  if !(new_mode in ("debug", "performance"))
+    throw(ArgumentError("Invalid execution mode: \"$(new_mode)\""))
+  end
+
+  # Set it in our runtime values, as well as saving it to disk
+  @set_preferences!("execution_mode" => new_mode)
+  @info("New execution mode set; restart your Julia session for this change to take effect!")
+end
+
+"""
+    set_debug_mode()
+
+  Equivalent to `set_execution_mode("debug")`.
+"""
+set_debug_mode() = set_execution_mode("debug")
+
+"""
+    set_performance_mode()
+
+  Equivalent to `set_execution_mode("performance")`.
+"""
+set_performance_mode() = set_execution_mode("performance")
+
+@static if VERSION >= v"1.6"
+  const execution_mode = @load_preference("execution_mode", "debug")
+else
+  const execution_mode = "debug"
+end

src/ODEs/ODESolvers/Tableaus.jl

@@ -221,11 +221,11 @@ struct IMEXTableau <: AbstractTableau{ImplicitExplicitTableau}
 
     msg = """Invalid IMEX tableau:
     the first tableau must be implicit and the second must be explicit."""
-    @assert (Tim <: ImplicitTableau && Tex == ExplicitTableau) msg
+    @check (Tim <: ImplicitTableau && Tex == ExplicitTableau) msg
 
     msg = """Invalid IMEX tableau:
     the nodes of the implicit and explicit tableaus must coincide."""
-    @assert isapprox(get_nodes(im_tableau), get_nodes(ex_tableau)) msg
+    @check isapprox(get_nodes(im_tableau), get_nodes(ex_tableau)) msg
 
     is_padded = _is_padded(im_tableau)
 

src/Polynomials/JacobiPolynomialBases.jl

@@ -47,7 +47,7 @@ return_type(::JacobiPolynomialBasis{D,T}) where {D,T} = T
 # Field implementation
 
 function return_cache(f::JacobiPolynomialBasis{D,T},x::AbstractVector{<:Point}) where {D,T}
-  @assert D == length(eltype(x)) "Incorrect number of point components"
+  @check D == length(eltype(x)) "Incorrect number of point components"
   np = length(x)
   ndof = length(f.terms)*num_components(T)
   n = 1 + _maximum(f.orders)

src/Polynomials/MonomialBases.jl

@@ -117,7 +117,7 @@ return_type(::MonomialBasis{D,T}) where {D,T} = T
 
 # Field implementation
 function return_cache(f::MonomialBasis{D,T},x::AbstractVector{<:Point}) where {D,T}
-  @assert D == length(eltype(x)) "Incorrect number of point components"
+  @check D == length(eltype(x)) "Incorrect number of point components"
   zT = zero(T)
   zxi = zero(eltype(eltype(x)))
   Tp = typeof( zT*zxi*zxi + zT*zxi*zxi  )
@@ -155,7 +155,7 @@ function _return_cache(
   TisbitsType::Val{true}) where {D,V,T}
 
   f = fg.fa
-  @assert D == length(eltype(x)) "Incorrect number of point components"
+  @check D == length(eltype(x)) "Incorrect number of point components"
   np = length(x)
   ndof = length(f.terms)*num_components(V)
   n = 1 + _maximum(f.orders)
@@ -253,7 +253,7 @@ function return_cache(
   x::AbstractVector{<:Point}) where {D,V}
 
   f = fg.fa
-  @assert D == length(eltype(x)) "Incorrect number of point components"
+  @check D == length(eltype(x)) "Incorrect number of point components"
   np = length(x)
   ndof = length(f.terms)*num_components(V)
   xi = testitem(x)

src/Polynomials/PCurlGradMonomialBases.jl

@@ -34,7 +34,7 @@ The `order` argument has the following meaning: the divergence of the  functions
 in this basis is in the P space of degree `order`.
 """
 function PCurlGradMonomialBasis{D}(::Type{T},order::Int) where {D,T}
-  @assert T<:Real "T needs to be <:Real since represents the type of the components of the vector value"
+  @check T<:Real "T needs to be <:Real since represents the type of the components of the vector value"
   P_k = MonomialBasis{D}(T, order, _p_filter)
   S_k = MonomialBasis{D}(T, order, _s_filter)
   pterms = P_k.terms
@@ -55,7 +55,7 @@ get_order(f::PCurlGradMonomialBasis{D,T}) where {D,T} = f.order
 return_type(::PCurlGradMonomialBasis{D,T}) where {D,T} = T
 
 function return_cache(f::PCurlGradMonomialBasis{D,T},x::AbstractVector{<:Point}) where {D,T}
-  @assert D == length(eltype(x)) "Incorrect number of point components"
+  @check D == length(eltype(x)) "Incorrect number of point components"
   np = length(x)
   ndof = _ndofs_pgrad(f)
   n = 1 + f.order+1
@@ -89,7 +89,7 @@ function return_cache(
   x::AbstractVector{<:Point})  where {D,T}
 
   f = fg.fa
-  @assert D == length(eltype(x)) "Incorrect number of point components"
+  @check D == length(eltype(x)) "Incorrect number of point components"
   np = length(x)
   ndof = _ndofs_pgrad(f)
   n = 1 + f.order+1

src/Polynomials/QCurlGradMonomialBases.jl

@@ -30,7 +30,7 @@ The `order` argument has the following meaning: the divergence of the  functions
 is in the Q space of degree `order`.
 """
 function QCurlGradMonomialBasis{D}(::Type{T},order::Int) where {D,T}
-  @assert T<:Real "T needs to be <:Real since represents the type of the components of the vector value"
+  @check T<:Real "T needs to be <:Real since represents the type of the components of the vector value"
   _order = order+1
   _t = tfill(_order,Val{D-1}())
   t = (_order+1,_t...)

src/Polynomials/QGradMonomialBases.jl

@@ -31,7 +31,7 @@ The `order` argument has the following meaning: the curl of the  functions in th
 is in the Q space of degree `order`.
 """
 function QGradMonomialBasis{D}(::Type{T},order::Int) where {D,T}
-  @assert T<:Real "T needs to be <:Real since represents the type of the components of the vector value"
+  @check T<:Real "T needs to be <:Real since represents the type of the components of the vector value"
   _order = order + 1
   _t = tfill(_order+1,Val{D-1}())
   t = (_order,_t...)
@@ -48,7 +48,7 @@ num_terms(f::QGradMonomialBasis{D,T}) where {D,T} = length(f.terms)*D
 get_order(f::QGradMonomialBasis) = f.order
 
 function return_cache(f::QGradMonomialBasis{D,T},x::AbstractVector{<:Point}) where {D,T}
-  @assert D == length(eltype(x)) "Incorrect number of point components"
+  @check D == length(eltype(x)) "Incorrect number of point components"
   np = length(x)
   ndof = _ndofs_qgrad(f)
   n = 1 + f.order+1
@@ -82,7 +82,7 @@ function return_cache(
   x::AbstractVector{<:Point}) where {D,T}
 
   f = fg.fa
-  @assert D == length(eltype(x)) "Incorrect number of point components"
+  @check D == length(eltype(x)) "Incorrect number of point components"
   np = length(x)
   ndof = _ndofs_qgrad(f)
   n = 1 + f.order+1

src/ReferenceFEs/BDMRefFEs.jl

@@ -122,7 +122,7 @@ function Conformity(reffe::GenericRefFE{BDM},sym::Symbol)
   function _BDM_face_values(p,et,order,phi)
 
     # Reference facet
-    @assert is_simplex(p) "We are assuming that all n-faces of the same n-dim are the same."
+    @check is_simplex(p) "We are assuming that all n-faces of the same n-dim are the same."
     fp = Polytope{num_dims(p)-1}(p,1)
 
     # geomap from ref face to polytope faces

src/ReferenceFEs/BezierRefFEs.jl

@@ -103,8 +103,8 @@ end
 ## Bernstein Basis
 
 function _bernstein_term(p,a,i)
-  @assert i ≤ p
-  @assert a ≤ p
+  @check i ≤ p
+  @check a ≤ p
   if ( i ≤ a ≤ p )
     f = factorial(p) ÷ ( factorial(i)*factorial(p-i) )
     b1 = binomial( p-i, a-i )
@@ -116,8 +116,8 @@ function _bernstein_term(p,a,i)
 end
 
 function _bernstein_term(p,a,b,i,j)
-  @assert i+j ≤ p
-  @assert a+b ≤ p
+  @check i+j ≤ p
+  @check a+b ≤ p
   if ( i ≤ a ≤ p-j ) && ( j ≤ b ≤ p-a )
     f = factorial(p) ÷ ( factorial(i)*factorial(j)*factorial(p-i-j) )
     b1 = binomial( p-i-j, a-i )
@@ -130,8 +130,8 @@ function _bernstein_term(p,a,b,i,j)
 end
 
 function _bernstein_term(p,a,b,c,i,j,k)
-  @assert i+j+k ≤ p
-  @assert a+b+c ≤ p
+  @check i+j+k ≤ p
+  @check a+b+c ≤ p
   if ( i ≤ a ≤ p-j-k ) && ( j ≤ b ≤ p-a-k ) && ( k ≤ c ≤ p-a-b )
     p!,i!,j!,k! = factorial(p),factorial(i),factorial(j),factorial(k)
     f = p! ÷ ( i!*j!*k!*factorial(p-i-j-k) )

src/ReferenceFEs/CDLagrangianRefFEs.jl

@@ -34,13 +34,13 @@ end
 
 function _CDLagrangianRefFE(::Type{T},p::ExtrusionPolytope{D},orders,cont) where {T,D}
   cond(c,o) = ( o > 0 || c == DISC )
-  @assert all((cond(cont[k],orders[k]) for k in 1:length(orders)))
+  @check all((cond(cont[k],orders[k]) for k in 1:length(orders)))
   _cd_lagrangian_ref_fe(T,p,orders,cont)
 end
 
 function _cd_lagrangian_ref_fe(::Type{T},p::ExtrusionPolytope{D},orders,cont) where {T,D}
 
-  @assert isa(p,ExtrusionPolytope)
+  @check isa(p,ExtrusionPolytope)
 
   prebasis = compute_monomial_basis(T,p,orders)