Untitled

using UnityEngine;

public static class Ballistics
{
    const float epsilon = 0.00000001f;

    public static Vector3 CalculateShootVelocity(Vector3 targetStartPos, Vector3 targetVelocity, Vector3 targetAcceleration, 
        Vector3 projectileStartPos, Vector3 projectileAcceleration, float projectileSpeed, out float duration)
    {
        Vector3 a = targetAcceleration - projectileAcceleration;
        Vector3 relativeTargetStartPos = targetStartPos - projectileStartPos;

        var roots = MathNet.Numerics.FindRoots.Polynomial(new double[]
        {
            Vector3.Dot(a, a) / 4.0f,
            Vector3.Dot(a, targetVelocity),
            Vector3.Dot(a, relativeTargetStartPos) + Vector3.SqrMagnitude(targetVelocity) - projectileSpeed * projectileSpeed,
            2.0f * Vector3.Dot(targetVelocity, relativeTargetStartPos),
            Vector3.SqrMagnitude(relativeTargetStartPos)
        });

        float? smallestPositiveRoot = null;
        foreach(var root in roots)
        {
            if (root.Imaginary != 0) continue;

            float rootReal = System.Convert.ToSingle(root.Real);
            if (rootReal > epsilon)
            {
                if (!smallestPositiveRoot.HasValue || rootReal < smallestPositiveRoot.Value) 
                    smallestPositiveRoot = rootReal;
            }
        }

        if(smallestPositiveRoot.HasValue)
        {
            duration = smallestPositiveRoot.Value;
            Vector3 travel = relativeTargetStartPos + targetVelocity * duration + duration * duration * (targetAcceleration / 2.0f);
            return travel / duration - (projectileAcceleration / 2.0f) * duration;
        }
        else
        {
            Debug.LogError("Unimplemented");
            duration = Mathf.Infinity;
            return Vector3.zero;
        }
    }
}