Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
14 changes: 13 additions & 1 deletion Sources/CelMek/NumericalSpacecraftPropagator.swift
Original file line number Diff line number Diff line change
Expand Up @@ -158,7 +158,19 @@ public struct NumericalSpacecraftPropagator: StateProvider {
let remaining = try forward ? target.duration(since: current) : current.duration(since: target)
let magnitude = min(nextStep,remaining.seconds)
let final = magnitude >= remaining.seconds
let end = final ? target : try current.advanced(by: Duration(seconds: forward ? magnitude : -magnitude))
var end = final ? target : try current.advanced(by: Duration(seconds: forward ? magnitude : -magnitude))
// Split epochs canonicalize their fractional seconds in 0..<1. A
// normal step can land a few rounding units before the endpoint
// (e.g. 0.1 + 0.5 on a 0.6 s query), leaving an unresolvable RK
// remainder. Integrate that rounding sliver with this step instead
// of dropping elapsed time or fabricating a state at a later epoch.
// Do not merge genuinely tiny requested/attempted steps: their stage
// resolution checks below must still reject collapsed timestamps.
let endpointRoundoffSeconds = 8*Double.ulpOfOne
if !final, magnitude > 16*endpointRoundoffSeconds,
abs(try target.duration(since:end).seconds) <= endpointRoundoffSeconds {
end = target
}
guard forward ? end > current : end < current else { throw NumericalPropagationError.epochResolutionExceeded }
let h = try end.duration(since: current).seconds
var stages: [[Double]] = [], candidate = y, depleted = false
Expand Down
18 changes: 18 additions & 0 deletions Tests/CelMekTests/NumericalPropagationTests.swift
Original file line number Diff line number Diff line change
Expand Up @@ -137,4 +137,22 @@ import Testing
}
}

@Test func fractionalStepRemaindersReachExactEndpointWithoutCollapsedStages() throws {
let frame=Frame(name:"fractional endpoint",pseudoInertial:true)
let dynamics=try SpacecraftDynamics(models:[Constant(acceleration:SIMD3(2,0,0))],context:ForceEvaluationContext(frame:frame,origin:.earth))
for fraction in [0.5,0.2693127683618036,0.99] {
for sign in [1.0,-1.0] {
let start=Epoch(offsetFromJ2000:try Duration(wholeSeconds:800000000,fractionalSecond:fraction))
let seed=SpacecraftState(cartesian:try CartesianState(epoch:start,frame:frame,position:.zero,velocity:SIMD3(3,0,0)),mass:try SpacecraftMassProperties(dryMassKilograms:1))
let options=try NumericalPropagationOptions(initialStepSeconds:0.1,maximumStepSeconds:2)
let model=try NumericalSpacecraftPropagator(initialState:seed,dynamics:dynamics,options:options)
let target=try start.advanced(by:Duration(seconds:sign*0.6))
let result=try model.propagation(at:target)
#expect(result.state.epoch==target)
#expect(abs(result.state.cartesian.position.x-(3*sign*0.6+0.36))<1e-12)
#expect(abs(result.state.cartesian.velocity.x-(3+2*sign*0.6))<1e-12)
}
}
}

}
Loading