Hi all,
I am working on fixing remaining bugs on the issue tracker prior to
9.0 release.
One of the issues is issue 324
<https://www.orekit.org/forge/issues/324> that
states that BoxAndSolarArraySpacecraft does not consider the lift
force. This
was first registered as a feature request, but I considered it a bug
and
changed its nature to bug some time ago.
Now it seems to me that the problem is *much* more serious than
initially
envisioned. I think BoxAndSolarArraySpacecraft not only forgot about
the
across-flux component, but it also computed the along-flux component
wrong!
I would like to have other developers have a look at my reasoning. It
is
really elementary physics, but I would like some confirmation, as the
current
code does not behave like this.
Let's consider a particle arriving at an angle on a plate and bouncing:
r
\ | plate
\ |
\|
n <--- /|
/ |
/ |
i |
The particle starts from the lower left, goes upward and to the right
along incident vector i, bounces, and changes its path upward and to
the
left along reflected vector r. The plate unit normal is n, here
oriented
to the left.
The effect of the drag on the plate is the opposite of the variation of
the linear momentum of the particle, and hence can be computed as
-mV (r - i). Considering an elastic collision, r - i is collinear to n.
If we consider the incident angle theta = pi - angle(i, n), we see
that r-i = 2cos(theta) n.
The force resulting from millions of particles hitting this plate will
be proportional to the cross section of the plate considering the flux
direction, which means this cross section will be A * cos(theta) is A
is the total surface of the plate.
The fact cos(theta) appears twice implies that the contribution of one
plate to the drag is proportional to the square of cos(theta), and is
always collinear to the plate normal.
The contribution of all the facets is only a summation of individual
contributions (we do not consider cast shadows in this class, this
would
imply building a full 3D model).
The previous simple elastic collision model is what I have implemented
to solve issue 324.
During validation of this new implementation, I was puzzled when
comparing with the previous implementation. After some head-scratching,
I now think this is because our former implementation was plain wrong.
Basically, what the former implementation did was compute the
cos(theta)
factor (using a simple dot product), and apply it to the incoming flux
in the flux direction. Even if we consider that the 2 in 2cos(theta)
can be removed from the model and included in the overall drag
coefficient
(which seems correct to me, so C will remain close to 2 which is often
a default value), there is a square missing and the direction is wrong.
Unfortunately, the wrong direction does not show up much in general
computation because when bodies are convex and symmetrical, plates on
the
left hand side and right hand side as seen from the incoming flux will
create forces cancelling out across-flux and adding together along
flux.
The global resultant is therefore really along flux (there is no lift
for symmetrical bodies). So we end up with the correct direction, but
is
is only by chance. However, the module of the contribution is wrong for
facets inclined. For head-on flux on a cubic body cos(theta)=1 and
there
is no error. For flux along the diagonal of a cubic body, the three
visible facets all have cos(theta)=1/sqrt(3), and the current code
uses a factor 0.577 when 0.333 should be used.
What do you think?
best regards,
Luc