From 81612761a374134a9c23175349d9019aa2138cbf Mon Sep 17 00:00:00 2001 From: AhmedSalih3d <36305327+AhmedSalih3d@users.noreply.github.com> Date: Fri, 5 Apr 2024 12:10:37 +0200 Subject: [PATCH] Functioning version --- example/StillWedge.jl | 406 ++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 406 insertions(+) create mode 100644 example/StillWedge.jl diff --git a/example/StillWedge.jl b/example/StillWedge.jl new file mode 100644 index 00000000..c405508d --- /dev/null +++ b/example/StillWedge.jl @@ -0,0 +1,406 @@ +using SPHExample +using StaticArrays +import StructArrays: StructArray +import LinearAlgebra: dot, norm, diagm, diag, cond, det +import Parameters: @unpack +import FastPow: @fastpow +import ProgressMeter: next!, finish! +using Format +using TimerOutputs +using Logging, LoggingExtras +using HDF5 +using Base.Threads + +function SaveVTKHDF(filepath,points, variable_names, args...) + @assert length(variable_names) == length(args) "Same number of variable_names as args is necessary" + h5open(filepath, "w") do io + # Create toplevel group /VTKHDF + gtop = HDF5.create_group(io, "VTKHDF") + + # Write version of VTKHDF format as an attribute + HDF5.attrs(gtop)["Version"] = [2, 1] + + # Write type of dataset ("PolyData") as an ASCII string to a "Type" attribute. + # This is a bit tricky because VTK/ParaView don't support UTF-8 strings here, which is + # the default in HDF5.jl. + let s = "PolyData" + dtype = HDF5.datatype(s) + HDF5.API.h5t_set_cset(dtype.id, HDF5.API.H5T_CSET_ASCII) + dspace = HDF5.dataspace(s) + attr = HDF5.create_attribute(gtop, "Type", dtype, dspace) + HDF5.write_attribute(attr, dtype, s) + end + + # Write points + number of points. + # Note that we need to reinterpret the vector of SVector onto a 3×Np matrix. + Np = length(points) + gtop["NumberOfPoints"] = [Np] + gtop["Points"] = reinterpret(reshape, eltype(eltype(points)), points) + + # Write velocities as point data. + let g = HDF5.create_group(gtop, "PointData") + for i ∈ eachindex(variable_names) + var_name = variable_names[i] + arg = args[i] + g[var_name] = reinterpret(reshape, eltype(eltype(arg)), arg) + end + end + + # Create and fill Vertices group. + let g = HDF5.create_group(gtop, "Vertices") + # In our case 1 point == 1 cell. + g["NumberOfCells"] = [Np] + g["NumberOfConnectivityIds"] = [Np] + g["Connectivity"] = collect(0:(Np - 1)) + g["Offsets"] = collect(0:Np) + close(g) + end + + # Add unused PolyData types. ParaView expects this, even if they're empty. + for type ∈ ("Lines", "Polygons", "Strips") + gempty = HDF5.create_group(gtop, type) + gempty["NumberOfCells"] = [0] + gempty["NumberOfConnectivityIds"] = [0] + gempty["Connectivity"] = Int[] + gempty["Offsets"] = [0] + close(gempty) + end + end +end + +# Really important to overload default function, gives 10x speed up? +# Overload the default function to do what you please +function ComputeInteractions!(SimMetaData, SimConstants, Position, KernelThreaded, KernelGradientThreaded, Density, Pressure, Velocity, dρdtI, dvdtI, i, j, MotionLimiter, ichunk) + @unpack FlagViscosityTreatment, FlagDensityDiffusion, FlagOutputKernelValues = SimMetaData + @unpack ρ₀, h, h⁻¹, m₀, αD, α, g, c₀, δᵩ, η², H², Cb⁻¹, ν₀, dx, SmagorinskyConstant, BlinConstant = SimConstants + + xᵢⱼ = Position[i] - Position[j] + xᵢⱼ² = dot(xᵢⱼ,xᵢⱼ) + if xᵢⱼ² <= H² + #https://discourse.julialang.org/t/sqrt-abs-x-is-even-faster-than-sqrt/58154/2 + dᵢⱼ = sqrt(abs(xᵢⱼ²)) + + q = min(dᵢⱼ * h⁻¹, 2.0) + invd²η² = inv(dᵢⱼ*dᵢⱼ+η²) + ∇ᵢWᵢⱼ = @fastpow (αD*5*(q-2)^3*q / (8h*(q*h+η²)) ) * xᵢⱼ + ρᵢ = Density[i] + ρⱼ = Density[j] + + vᵢ = Velocity[i] + vⱼ = Velocity[j] + vᵢⱼ = vᵢ - vⱼ + density_symmetric_term = dot(-vᵢⱼ, ∇ᵢWᵢⱼ) + dρdt⁺ = - ρᵢ * (m₀/ρⱼ) * density_symmetric_term + dρdt⁻ = - ρⱼ * (m₀/ρᵢ) * density_symmetric_term + + # Density diffusion + if FlagDensityDiffusion + Pᵢⱼᴴ = ρ₀ * (-g) * -xᵢⱼ[end] + ρᵢⱼᴴ = InverseHydrostaticEquationOfState(ρ₀, Pᵢⱼᴴ, Cb⁻¹) + Pⱼᵢᴴ = -Pᵢⱼᴴ + ρⱼᵢᴴ = InverseHydrostaticEquationOfState(ρ₀, Pⱼᵢᴴ, Cb⁻¹) + + ρⱼᵢ = ρⱼ - ρᵢ + MLcond = MotionLimiter[i] * MotionLimiter[j] + ddt_symmetric_term = δᵩ * h * c₀ * 2 * invd²η² * dot(-xᵢⱼ, ∇ᵢWᵢⱼ) #* MLcond + Dᵢ = ddt_symmetric_term * (m₀/ρⱼ) * ( ρⱼᵢ - ρᵢⱼᴴ) * MLcond + Dⱼ = ddt_symmetric_term * (m₀/ρᵢ) * (-ρⱼᵢ - ρⱼᵢᴴ) * MLcond + else + Dᵢ = 0.0 + Dⱼ = 0.0 + end + dρdtI[ichunk][i] += dρdt⁺ + Dᵢ + dρdtI[ichunk][j] += dρdt⁻ + Dⱼ + + + Pᵢ = Pressure[i] + Pⱼ = Pressure[j] + Pfac = (Pᵢ+Pⱼ)/(ρᵢ*ρⱼ) + dvdt⁺ = - m₀ * Pfac * ∇ᵢWᵢⱼ + dvdt⁻ = - dvdt⁺ + + if FlagViscosityTreatment == :ArtificialViscosity + ρ̄ᵢⱼ = (ρᵢ+ρⱼ)*0.5 + cond = dot(vᵢⱼ, xᵢⱼ) + cond_bool = cond < 0.0 + μᵢⱼ = h*cond * invd²η² + Πᵢ = - m₀ * (cond_bool*(-α*c₀*μᵢⱼ)/ρ̄ᵢⱼ) * ∇ᵢWᵢⱼ + Πⱼ = - Πᵢ + else + Πᵢ = zero(xᵢⱼ) + Πⱼ = Πᵢ + end + + if FlagViscosityTreatment == :Laminar || FlagViscosityTreatment == :LaminarSPS + # 4 comes from 2 divided by 0.5 from average density + # should divide by ρᵢ eq 6 DPC + # ν₀∇²uᵢ = (1/ρᵢ) * ( (4 * m₀ * (ρᵢ * ν₀) * dot( xᵢⱼ, ∇ᵢWᵢⱼ) ) / ( (ρᵢ + ρⱼ) + (dᵢⱼ * dᵢⱼ + η²) ) ) * vᵢⱼ + # ν₀∇²uⱼ = (1/ρⱼ) * ( (4 * m₀ * (ρⱼ * ν₀) * dot(-xᵢⱼ,-∇ᵢWᵢⱼ) ) / ( (ρᵢ + ρⱼ) + (dᵢⱼ * dᵢⱼ + η²) ) ) * -vᵢⱼ + visc_symmetric_term = (4 * m₀ * ν₀ * dot( xᵢⱼ, ∇ᵢWᵢⱼ)) / ((ρᵢ + ρⱼ) + (dᵢⱼ * dᵢⱼ + η²)) + # ν₀∇²uᵢ = (1/ρᵢ) * visc_symmetric_term * vᵢⱼ * ρᵢ + # ν₀∇²uⱼ = (1/ρⱼ) * visc_symmetric_term * -vᵢⱼ * ρⱼ + ν₀∇²uᵢ = visc_symmetric_term * vᵢⱼ + ν₀∇²uⱼ = -ν₀∇²uᵢ #visc_symmetric_term * -vᵢⱼ + else + ν₀∇²uᵢ = zero(xᵢⱼ) + ν₀∇²uⱼ = ν₀∇²uᵢ + end + + if FlagViscosityTreatment == :LaminarSPS + Iᴹ = diagm(one.(xᵢⱼ)) + #julia> a .- a' + # 3×3 SMatrix{3, 3, Float64, 9} with indices SOneTo(3)×SOneTo(3): + # 0.0 0.0 0.0 + # 0.0 0.0 0.0 + # 0.0 0.0 0.0 + # Strain *rate* tensor is the gradient of velocity + Sᵢ = ∇vᵢ = (m₀/ρⱼ) * (vⱼ - vᵢ) * ∇ᵢWᵢⱼ' + norm_Sᵢ = sqrt(2 * sum(Sᵢ .^ 2)) + νtᵢ = (SmagorinskyConstant * dx)^2 * norm_Sᵢ + trace_Sᵢ = sum(diag(Sᵢ)) + τᶿᵢ = 2*νtᵢ*ρᵢ * (Sᵢ - (1/3) * trace_Sᵢ * Iᴹ) - (2/3) * ρᵢ * BlinConstant * dx^2 * norm_Sᵢ^2 * Iᴹ + Sⱼ = ∇vⱼ = (m₀/ρᵢ) * (vᵢ - vⱼ) * -∇ᵢWᵢⱼ' + norm_Sⱼ = sqrt(2 * sum(Sⱼ .^ 2)) + νtⱼ = (SmagorinskyConstant * dx)^2 * norm_Sⱼ + trace_Sⱼ = sum(diag(Sⱼ)) + τᶿⱼ = 2*νtⱼ*ρⱼ * (Sⱼ - (1/3) * trace_Sⱼ * Iᴹ) - (2/3) * ρⱼ * BlinConstant * dx^2 * norm_Sⱼ^2 * Iᴹ + + # MATHEMATICALLY THIS IS DOT PRODUCT TO GO FROM TENSOR TO VECTOR, BUT USE * IN JULIA TO REPRESENT IT + dτdtᵢ = (m₀/(ρⱼ * ρᵢ)) * (τᶿᵢ + τᶿⱼ) * ∇ᵢWᵢⱼ + dτdtⱼ = (m₀/(ρᵢ * ρⱼ)) * (τᶿᵢ + τᶿⱼ) * -∇ᵢWᵢⱼ + else + dτdtᵢ = zero(xᵢⱼ) + dτdtⱼ = dτdtᵢ + end + + dvdtI[ichunk][i] += dvdt⁺ + Πᵢ + ν₀∇²uᵢ + dτdtᵢ + dvdtI[ichunk][j] += dvdt⁻ + Πⱼ + ν₀∇²uⱼ + dτdtⱼ + + if FlagOutputKernelValues + Wᵢⱼ = @fastpow αD*(1-q/2)^4*(2*q + 1) + KernelThreaded[ichunk][i] += Wᵢⱼ + KernelThreaded[ichunk][j] += Wᵢⱼ + KernelGradientThreaded[ichunk][i] += ∇ᵢWᵢⱼ + KernelGradientThreaded[ichunk][j] += -∇ᵢWᵢⱼ + end + end + + return nothing +end + +# Reduce threaded arrays +@inline function reduce_sum!(target_array, arrays) + for array in arrays + target_array .+= array + end +end +@inbounds function SimulationLoop(ComputeInteractions!, SimMetaData, SimConstants, SimParticles, Stencil, ParticleRanges, UniqueCells, SortingScratchSpace, Kernel, KernelThreaded, KernelGradient, KernelGradientThreaded, dρdtI, dρdtIThreaded, AccelerationThreaded, Velocityₙ⁺, Positionₙ⁺, ρₙ⁺, InverseCutOff) + Position = SimParticles.Position + Density = SimParticles.Density + Pressure = SimParticles.Pressure + Velocity = SimParticles.Velocity + Acceleration = SimParticles.Acceleration + GravityFactor = SimParticles.GravityFactor + MotionLimiter = SimParticles.MotionLimiter + + @timeit SimMetaData.HourGlass "01 Update TimeStep" dt = Δt(Position, Velocity, Acceleration, SimConstants) + dt₂ = dt * 0.5 + + # if mod(SimMetaData.Iteration,10) == 0 + @timeit SimMetaData.HourGlass "02 Calculate IndexCounter" IndexCounter = UpdateNeighbors!(SimParticles, InverseCutOff, SortingScratchSpace, ParticleRanges, UniqueCells) + # else + # IndexCounter = findfirst(isequal(0), ParticleRanges) - 2 + # end + + @timeit SimMetaData.HourGlass "03 ResetArrays" ResetArrays!(Kernel, KernelGradient, dρdtI, Acceleration); ResetArrays!.(KernelThreaded, KernelGradientThreaded, dρdtIThreaded, AccelerationThreaded) + + Pressure!(SimParticles.Pressure,SimParticles.Density,SimConstants) + @timeit SimMetaData.HourGlass "04 First NeighborLoop" NeighborLoop!(ComputeInteractions!, SimMetaData, SimConstants, ParticleRanges, Stencil, Position, KernelThreaded, KernelGradientThreaded, Density, Pressure, Velocity, dρdtIThreaded, AccelerationThreaded, MotionLimiter, UniqueCells, IndexCounter) + @timeit SimMetaData.HourGlass "04A Reduction" reduce_sum!(dρdtI, dρdtIThreaded) + @timeit SimMetaData.HourGlass "04B Reduction" reduce_sum!(Acceleration, AccelerationThreaded) + + @timeit SimMetaData.HourGlass "05 Update To Half TimeStep" @inbounds for i in eachindex(Position) + Acceleration[i] += ConstructGravitySVector(Acceleration[i], SimConstants.g * GravityFactor[i]) + Positionₙ⁺[i] = Position[i] + Velocity[i] * dt₂ * MotionLimiter[i] + Velocityₙ⁺[i] = Velocity[i] + Acceleration[i] * dt₂ * MotionLimiter[i] + ρₙ⁺[i] = Density[i] + dρdtI[i] * dt₂ + end + + @timeit SimMetaData.HourGlass "06 Half LimitDensityAtBoundary" LimitDensityAtBoundary!(ρₙ⁺, SimConstants.ρ₀, MotionLimiter) + + @timeit SimMetaData.HourGlass "07 ResetArrays" ResetArrays!(Kernel, KernelGradient, dρdtI, Acceleration); ResetArrays!.(KernelThreaded, KernelGradientThreaded, dρdtIThreaded, AccelerationThreaded) + + Pressure!(SimParticles.Pressure, ρₙ⁺,SimConstants) + @timeit SimMetaData.HourGlass "08 Second NeighborLoop" NeighborLoop!(ComputeInteractions!, SimMetaData, SimConstants, ParticleRanges, Stencil, Positionₙ⁺, KernelThreaded, KernelGradientThreaded, ρₙ⁺, Pressure, Velocityₙ⁺, dρdtIThreaded, AccelerationThreaded, MotionLimiter, UniqueCells, IndexCounter) + @timeit SimMetaData.HourGlass "08A Reduction" reduce_sum!(dρdtI, dρdtIThreaded) + @timeit SimMetaData.HourGlass "08B Reduction" reduce_sum!(Acceleration, AccelerationThreaded) + + @timeit SimMetaData.HourGlass "09 Final LimitDensityAtBoundary" LimitDensityAtBoundary!(Density, SimConstants.ρ₀, MotionLimiter) + + @timeit SimMetaData.HourGlass "10 Final Density" DensityEpsi!(Density, dρdtI, ρₙ⁺, dt) + + + @timeit SimMetaData.HourGlass "11 Update To Final TimeStep" @inbounds for i in eachindex(Position) + Acceleration[i] += ConstructGravitySVector(Acceleration[i], SimConstants.g * GravityFactor[i]) + Velocity[i] += Acceleration[i] * dt * MotionLimiter[i] + Position[i] += (((Velocity[i] + (Velocity[i] - Acceleration[i] * dt * MotionLimiter[i])) / 2) * dt) * MotionLimiter[i] + end + + SimMetaData.Iteration += 1 + SimMetaData.CurrentTimeStep = dt + SimMetaData.TotalTime += dt + + if SimMetaData.FlagOutputKernelValues + reduce_sum!(Kernel, KernelThreaded) + reduce_sum!(KernelGradient, KernelGradientThreaded) + end + + return nothing +end + +###=== +function RunSimulation(;FluidCSV::String, + BoundCSV::String, + SimMetaData::SimulationMetaData{Dimensions, FloatType}, + SimConstants::SimulationConstants, + SimLogger::SimulationLogger + ) where {Dimensions,FloatType} + + if SimMetaData.FlagLog + InitializeLogger(SimLogger,SimConstants,SimMetaData) + end + + # If save directory is not already made, make it + if !isdir(SimMetaData.SaveLocation) + mkdir(SimMetaData.SaveLocation) + end + + # Delete previous result files + foreach(rm, filter(endswith(".vtp"), readdir(SimMetaData.SaveLocation,join=true))) + try + foreach(rm, filter(endswith(".h5"), readdir(SimMetaData.SaveLocation,join=true))) + catch + end + + # Unpack the relevant simulation meta data + @unpack HourGlass, SaveLocation, SimulationName, SilentOutput, ThreadsCPU = SimMetaData; + + # Load in particles + SimParticles, dρdtI, Velocityₙ⁺, Positionₙ⁺, ρₙ⁺, Kernel, KernelGradient = AllocateDataStructures(Dimensions,FloatType, FluidCSV,BoundCSV) + + + NumberOfPoints = length(SimParticles)::Int #Have to type declare, else error? + @inline Pressure!(SimParticles.Pressure,SimParticles.Density,SimConstants) + + @inline begin + n_copy = Base.Threads.nthreads() + KernelThreaded = [copy(Kernel) for _ in 1:n_copy] + KernelGradientThreaded = [copy(KernelGradient) for _ in 1:n_copy] + dρdtIThreaded = [copy(dρdtI) for _ in 1:n_copy] + AccelerationThreaded = [copy(KernelGradient) for _ in 1:n_copy] + end + + # Produce sorting related variables + ParticleRanges = zeros(Int, NumberOfPoints + 1) + UniqueCells = zeros(CartesianIndex{Dimensions}, NumberOfPoints) + Stencil = ConstructStencil(Val(Dimensions)) + _, SortingScratchSpace = Base.Sort.make_scratch(nothing, eltype(SimParticles), NumberOfPoints) + + # Produce data saving functions + SaveLocation_ = SimMetaData.SaveLocation * "/" * SimulationName * "_" * lpad(SimMetaData.Iteration,6,"0") * ".vtkhdf" + # FidBig = h5open(SaveLocation_, "w") + + # SaveFile = (GroupName) -> SaveHDF5!(FidBig, GroupName, ["Position", "Kernel", "KernelGradient", "Density", "Pressure","Velocity", "Acceleration", "BoundaryBool" , "ID"], to_3d(SimParticles.Position), Kernel, KernelGradient, SimParticles.Density, SimParticles.Pressure, SimParticles.Velocity, SimParticles.Acceleration, Int.(SimParticles.BoundaryBool), SimParticles.ID) + # SaveFile = (GroupName) -> SaveHDF5OneBig!(FidBig, GroupName, ["Position", "Kernel", "KernelGradient", "Density", "Pressure","Velocity", "Acceleration", "BoundaryBool" , "ID"], to_3d(SimParticles.Position), Kernel, KernelGradient, SimParticles.Density, SimParticles.Pressure, SimParticles.Velocity, SimParticles.Acceleration, SimParticles.ID) + SaveFile = (SaveLocation_) -> SaveVTKHDF(SaveLocation_,to_3d(SimParticles.Position),["Kernel", "KernelGradient", "Density", "Pressure","Velocity", "Acceleration", "BoundaryBool" , "ID"], Kernel, KernelGradient, SimParticles.Density, SimParticles.Pressure, SimParticles.Velocity, SimParticles.Acceleration, Int.(SimParticles.BoundaryBool), SimParticles.ID) + @inline SaveFile(SaveLocation_) + SimMetaData.OutputIterationCounter += 1 #Since a file has been saved + + + InverseCutOff = Val(1/(SimConstants.H)) + + # Normal run and save data + generate_showvalues(Iteration, TotalTime, TimeLeftInSeconds) = () -> [(:(Iteration),format(FormatExpr("{1:d}"), Iteration)), (:(TotalTime),format(FormatExpr("{1:3.3f}"), TotalTime)), (:(TimeLeftInSeconds),format(FormatExpr("{1:3.1f} [s]"), TimeLeftInSeconds))] + + # This is for some reason to trick the compiler to avoid dispatch error on SimulationLoop due to SimParticles + # @inline on SimulationLoop directly slows down code + f = () -> SimulationLoop(ComputeInteractions!, SimMetaData, SimConstants, SimParticles, Stencil, ParticleRanges, UniqueCells, SortingScratchSpace, Kernel, KernelThreaded, KernelGradient, KernelGradientThreaded, dρdtI, dρdtIThreaded, AccelerationThreaded, Velocityₙ⁺, Positionₙ⁺, ρₙ⁺, InverseCutOff) + + @inbounds while true + + @inline f() + + if SimMetaData.TotalTime >= SimMetaData.OutputEach * SimMetaData.OutputIterationCounter + + # @timeit HourGlass "12A Output Data" SaveFile(string(SimMetaData.OutputIterationCounter)) + SaveLocation_ = SimMetaData.SaveLocation * "/" * SimulationName * "_" * lpad(SimMetaData.Iteration,6,"0") * ".vtkhdf" + @timeit HourGlass "12A Output Data" SaveFile(SaveLocation_) + + if SimMetaData.FlagLog + LogStep(SimLogger, SimMetaData, HourGlass) + SimMetaData.StepsTakenForLastOutput = SimMetaData.Iteration + end + + SimMetaData.OutputIterationCounter += 1 + end + + if !SilentOutput + TimeLeftInSeconds = (SimMetaData.SimulationTime - SimMetaData.TotalTime) * (TimerOutputs.tottime(HourGlass)/1e9 / SimMetaData.TotalTime) + @timeit HourGlass "13 Next TimeStep" next!(SimMetaData.ProgressSpecification; showvalues = generate_showvalues(SimMetaData.Iteration , SimMetaData.TotalTime, TimeLeftInSeconds)) + end + + if SimMetaData.TotalTime > SimMetaData.SimulationTime + # @timeit HourGlass "12B Close h5 output file" close(FidBig) + + finish!(SimMetaData.ProgressSpecification) + show(HourGlass,sortby=:name) + show(HourGlass) + + if SimMetaData.FlagLog + LogFinal(SimLogger, HourGlass) + close(SimLogger.LoggerIo) + end + + break + end + end +end + + +let + Dimensions = 2 + FloatType = Float64 + + SimMetaDataWedge = SimulationMetaData{Dimensions,FloatType}( + SimulationName="Test", + SaveLocation="E:/SecondApproach/TESTING_CPU", + SimulationTime=4, + OutputEach=0.01, + FlagDensityDiffusion=true, + FlagOutputKernelValues=false, + FlagLog=true + ) + + SimConstantsWedge = SimulationConstants{FloatType}(dx=0.02,c₀=42.48576250492629, δᵩ = 0.1, CFL=0.2) + + SimLogger = SimulationLogger(SimMetaDataWedge.SaveLocation) + + # Remove '@profview' if you do not want VS Code timers + # println(@report_call target_modules=(@__MODULE__,) RunSimulation( + # FluidCSV = "./input/still_wedge/StillWedge_Dp0.02_Fluid.csv", + # BoundCSV = "./input/still_wedge/StillWedge_Dp0.02_Bound.csv", + # SimMetaData = SimMetaDataWedge, + # SimConstants = SimConstantsWedge, + # SimLogger = SimLogger, + # )) + + RunSimulation( + FluidCSV = "./input/still_wedge/StillWedge_Dp0.02_Fluid.csv", + BoundCSV = "./input/still_wedge/StillWedge_Dp0.02_Bound.csv", + # FluidCSV = "./input/dam_break_2d/DamBreak2d_Dp0.02_Fluid.csv", + # BoundCSV = "./input/dam_break_2d/DamBreak2d_Dp0.02_Bound.csv", + SimMetaData = SimMetaDataWedge, + SimConstants = SimConstantsWedge, + SimLogger = SimLogger + ) +end