Mercurial > repos > public > sbplib_julia
comparison test/Grids/mapped_grid_test.jl @ 1777:819ab806960f feature/grids/manifolds
Merge feature/grids/curvilinear
| author | Jonatan Werpers <jonatan@werpers.com> |
|---|---|
| date | Sun, 15 Sep 2024 23:00:15 +0200 |
| parents | 03894fd7b132 ecec2b0eea0f |
| children | 2b5f81e288f1 |
comparison
equal
deleted
inserted
replaced
| 1750:e2abd72d7ce0 | 1777:819ab806960f |
|---|---|
| 35 J = map(ξ̄ -> @SArray(fill(2., 2, 2)), lg) | 35 J = map(ξ̄ -> @SArray(fill(2., 2, 2)), lg) |
| 36 mg = MappedGrid(lg, x̄, J) | 36 mg = MappedGrid(lg, x̄, J) |
| 37 | 37 |
| 38 @test mg isa Grid{SVector{2, Float64},2} | 38 @test mg isa Grid{SVector{2, Float64},2} |
| 39 @test jacobian(mg) isa Array{<:AbstractMatrix} | 39 @test jacobian(mg) isa Array{<:AbstractMatrix} |
| 40 @test logicalgrid(mg) isa Grid | 40 @test logical_grid(mg) isa Grid |
| 41 | 41 |
| 42 @test collect(mg) == x̄ | 42 @test collect(mg) == x̄ |
| 43 @test jacobian(mg) == J | 43 @test jacobian(mg) == J |
| 44 @test logicalgrid(mg) == lg | 44 @test logical_grid(mg) == lg |
| 45 | 45 |
| 46 | 46 |
| 47 x̄ = map(ξ̄ -> @SVector[ξ̄[1],ξ̄[2], ξ̄[1] + ξ̄[2]], lg) | 47 x̄ = map(ξ̄ -> @SVector[ξ̄[1],ξ̄[2], ξ̄[1] + ξ̄[2]], lg) |
| 48 J = map(ξ̄ -> @SMatrix[1 0; 0 1; 1 1], lg) | 48 J = map(ξ̄ -> @SMatrix[1 0; 0 1; 1 1], lg) |
| 49 mg = MappedGrid(lg, x̄, J) | 49 mg = MappedGrid(lg, x̄, J) |
| 50 | 50 |
| 51 @test mg isa Grid{SVector{3, Float64},2} | 51 @test mg isa Grid{SVector{3, Float64},2} |
| 52 @test jacobian(mg) isa Array{<:AbstractMatrix} | 52 @test jacobian(mg) isa Array{<:AbstractMatrix} |
| 53 @test logicalgrid(mg) isa Grid | 53 @test logical_grid(mg) isa Grid |
| 54 | 54 |
| 55 @test collect(mg) == x̄ | 55 @test collect(mg) == x̄ |
| 56 @test jacobian(mg) == J | 56 @test jacobian(mg) == J |
| 57 @test logicalgrid(mg) == lg | 57 @test logical_grid(mg) == lg |
| 58 | 58 |
| 59 sz1 = (10,11) | 59 sz1 = (10,11) |
| 60 sz2 = (10,12) | 60 sz2 = (10,12) |
| 61 @test_throws ArgumentError("Sizes must match") MappedGrid( | 61 @test_throws ArgumentError("Sizes must match") MappedGrid( |
| 62 equidistant_grid((0,0), (1,1), sz2...), | 62 equidistant_grid((0,0), (1,1), sz2...), |
| 245 0; | 245 0; |
| 246 1+ξ*(ξ-1); | 246 1+ξ*(ξ-1); |
| 247 ] | 247 ] |
| 248 | 248 |
| 249 function expected_bg(mg, bId, Jb) | 249 function expected_bg(mg, bId, Jb) |
| 250 lg = logicalgrid(mg) | 250 lg = logical_grid(mg) |
| 251 return MappedGrid( | 251 return MappedGrid( |
| 252 boundary_grid(lg, bId), | 252 boundary_grid(lg, bId), |
| 253 map(x̄, boundary_grid(lg, bId)), | 253 map(x̄, boundary_grid(lg, bId)), |
| 254 map(Jb, boundary_grid(lg, bId)), | 254 map(Jb, boundary_grid(lg, bId)), |
| 255 ) | 255 ) |
| 277 x̄, J = _partially_curved_mapping() | 277 x̄, J = _partially_curved_mapping() |
| 278 mg = mapped_grid(x̄, J, 10, 11) | 278 mg = mapped_grid(x̄, J, 10, 11) |
| 279 @test mg isa MappedGrid{SVector{2,Float64}, 2} | 279 @test mg isa MappedGrid{SVector{2,Float64}, 2} |
| 280 | 280 |
| 281 lg = equidistant_grid((0,0), (1,1), 10, 11) | 281 lg = equidistant_grid((0,0), (1,1), 10, 11) |
| 282 @test logicalgrid(mg) == lg | 282 @test logical_grid(mg) == lg |
| 283 @test collect(mg) == map(x̄, lg) | 283 @test collect(mg) == map(x̄, lg) |
| 284 | 284 |
| 285 @test mapped_grid(lg, x̄, J) == mg | 285 @test mapped_grid(lg, x̄, J) == mg |
| 286 end | |
| 287 | |
| 288 @testset "jacobian_determinant" begin | |
| 289 x̄((ξ, η)) = @SVector[ξ*η, ξ + η^2] | |
| 290 J((ξ, η)) = @SMatrix[ | |
| 291 η ξ; | |
| 292 1 2η; | |
| 293 ] | |
| 294 | |
| 295 g = mapped_grid(x̄, J, 10, 11) | |
| 296 J = map(logicalgrid(g)) do (ξ,η) | |
| 297 2η^2 - ξ | |
| 298 end | |
| 299 @test jacobian_determinant(g) ≈ J | |
| 300 | |
| 301 | |
| 302 lg = equidistant_grid((0,0), (1,1), 11, 21) | |
| 303 x̄ = map(ξ̄ -> @SVector[ξ̄[1],ξ̄[2], ξ̄[1] + ξ̄[2]], lg) | |
| 304 J = map(ξ̄ -> @SMatrix[1 0; 0 1; 1 1], lg) | |
| 305 mg = MappedGrid(lg, x̄, J) | |
| 306 | |
| 307 @test_broken jacobian(mg) isa AbstractArray{2,Float64} | |
| 308 end | 286 end |
| 309 | 287 |
| 310 @testset "metric_tensor" begin | 288 @testset "metric_tensor" begin |
| 311 x̄((ξ, η)) = @SVector[ξ*η, ξ + η^2] | 289 x̄((ξ, η)) = @SVector[ξ*η, ξ + η^2] |
| 312 J((ξ, η)) = @SMatrix[ | 290 J((ξ, η)) = @SMatrix[ |
| 313 η ξ; | 291 η ξ; |
| 314 1 2η; | 292 1 2η; |
| 315 ] | 293 ] |
| 316 | 294 |
| 317 g = mapped_grid(x̄, J, 10, 11) | 295 g = mapped_grid(x̄, J, 10, 11) |
| 318 G = map(logicalgrid(g)) do (ξ,η) | 296 G = map(logical_grid(g)) do (ξ,η) |
| 319 @SMatrix[ | 297 @SMatrix[ |
| 320 1+η^2 ξ*η+2η; | 298 1+η^2 ξ*η+2η; |
| 321 ξ*η+2η ξ^2 + 4η^2; | 299 ξ*η+2η ξ^2 + 4η^2; |
| 322 ] | 300 ] |
| 323 end | 301 end |
| 330 1+ξ 0; | 308 1+ξ 0; |
| 331 ξ 1+η; | 309 ξ 1+η; |
| 332 ] | 310 ] |
| 333 | 311 |
| 334 g = mapped_grid(x̄, J, 10, 11) | 312 g = mapped_grid(x̄, J, 10, 11) |
| 335 G⁻¹ = map(logicalgrid(g)) do (ξ,η) | 313 G⁻¹ = map(logical_grid(g)) do (ξ,η) |
| 336 @SMatrix[ | 314 @SMatrix[ |
| 337 (1+η)^2 -ξ*(1+η); | 315 (1+η)^2 -ξ*(1+η); |
| 338 -ξ*(1+η) (1+ξ)^2+ξ^2; | 316 -ξ*(1+η) (1+ξ)^2+ξ^2; |
| 339 ]/(((1+ξ)^2+ξ^2)*(1+η)^2 - ξ^2*(1+η)^2) | 317 ]/(((1+ξ)^2+ξ^2)*(1+η)^2 - ξ^2*(1+η)^2) |
| 340 | 318 |
| 382 g = mapped_grid(_partially_curved_mapping()...,10, 11) | 360 g = mapped_grid(_partially_curved_mapping()...,10, 11) |
| 383 | 361 |
| 384 @test normal(g, CartesianBoundary{1,LowerBoundary}()) == fill(@SVector[-1,0], 11) | 362 @test normal(g, CartesianBoundary{1,LowerBoundary}()) == fill(@SVector[-1,0], 11) |
| 385 @test normal(g, CartesianBoundary{1,UpperBoundary}()) == fill(@SVector[1,0], 11) | 363 @test normal(g, CartesianBoundary{1,UpperBoundary}()) == fill(@SVector[1,0], 11) |
| 386 @test normal(g, CartesianBoundary{2,LowerBoundary}()) == fill(@SVector[0,-1], 10) | 364 @test normal(g, CartesianBoundary{2,LowerBoundary}()) == fill(@SVector[0,-1], 10) |
| 387 @test normal(g, CartesianBoundary{2,UpperBoundary}()) ≈ map(boundary_grid(g,CartesianBoundary{2,UpperBoundary}())|>logicalgrid) do ξ̄ | 365 @test normal(g, CartesianBoundary{2,UpperBoundary}()) ≈ map(boundary_grid(g,CartesianBoundary{2,UpperBoundary}())|>logical_grid) do ξ̄ |
| 388 α = 1-2ξ̄[1] | 366 α = 1-2ξ̄[1] |
| 389 @SVector[α,1]/√(α^2 + 1) | 367 @SVector[α,1]/√(α^2 + 1) |
| 390 end | 368 end |
| 391 | 369 |
| 392 g = mapped_grid(_fully_curved_mapping()...,5,4) | 370 g = mapped_grid(_fully_curved_mapping()...,5,4) |
| 393 | 371 |
| 394 unit(v) = v/norm(v) | 372 unit(v) = v/norm(v) |
| 395 @testset let bId = CartesianBoundary{1,LowerBoundary}() | 373 @testset let bId = CartesianBoundary{1,LowerBoundary}() |
| 396 lbg = boundary_grid(logicalgrid(g), bId) | 374 lbg = boundary_grid(logical_grid(g), bId) |
| 397 @test normal(g, bId) ≈ map(lbg) do (ξ, η) | 375 @test normal(g, bId) ≈ map(lbg) do (ξ, η) |
| 398 -unit(@SVector[1/2, η/3-1/6]) | 376 -unit(@SVector[1/2, η/3-1/6]) |
| 399 end | 377 end |
| 400 end | 378 end |
| 401 | 379 |
| 402 @testset let bId = CartesianBoundary{1,UpperBoundary}() | 380 @testset let bId = CartesianBoundary{1,UpperBoundary}() |
| 403 lbg = boundary_grid(logicalgrid(g), bId) | 381 lbg = boundary_grid(logical_grid(g), bId) |
| 404 @test normal(g, bId) ≈ map(lbg) do (ξ, η) | 382 @test normal(g, bId) ≈ map(lbg) do (ξ, η) |
| 405 unit(@SVector[7/2, 2η-1]/(5 + 3η + 2η^2)) | 383 unit(@SVector[7/2, 2η-1]/(5 + 3η + 2η^2)) |
| 406 end | 384 end |
| 407 end | 385 end |
| 408 | 386 |
| 409 @testset let bId = CartesianBoundary{2,LowerBoundary}() | 387 @testset let bId = CartesianBoundary{2,LowerBoundary}() |
| 410 lbg = boundary_grid(logicalgrid(g), bId) | 388 lbg = boundary_grid(logical_grid(g), bId) |
| 411 @test normal(g, bId) ≈ map(lbg) do (ξ, η) | 389 @test normal(g, bId) ≈ map(lbg) do (ξ, η) |
| 412 -unit(@SVector[-2ξ, 2]/(6 + ξ^2 - 2ξ)) | 390 -unit(@SVector[-2ξ, 2]/(6 + ξ^2 - 2ξ)) |
| 413 end | 391 end |
| 414 end | 392 end |
| 415 | 393 |
| 416 @testset let bId = CartesianBoundary{2,UpperBoundary}() | 394 @testset let bId = CartesianBoundary{2,UpperBoundary}() |
| 417 lbg = boundary_grid(logicalgrid(g), bId) | 395 lbg = boundary_grid(logical_grid(g), bId) |
| 418 @test normal(g, bId) ≈ map(lbg) do (ξ, η) | 396 @test normal(g, bId) ≈ map(lbg) do (ξ, η) |
| 419 unit(@SVector[-3ξ, 2]/(6 + ξ^2 + 3ξ)) | 397 unit(@SVector[-3ξ, 2]/(6 + ξ^2 + 3ξ)) |
| 420 end | 398 end |
| 421 end | 399 end |
| 422 end | 400 end |
