Diophantine linear equation solver in C#

Diophantine linear equations are equations of two unkowns on the format

ax + by = c

The code below shows how we can solve Diophantine equations with a sample equation :

7x - 9y = 3

DiophantineSolver.cs

class Program
{
	static void Main()
	{
		int a = 7,
		b = -9, 
		c = 3;
		DiophantineSolver.DiophantineSolution solution = DiophantineSolver.ExtendedGcd(a, b, c);

		Console.WriteLine(solution.ToString());
	}
}

/// <summary>
/// Class containing methods to solve Diophantine equations using the Extended Euclidean Algorithm.
/// </summary>
class DiophantineSolver
{
	/// <summary>
	/// Struct to hold the solution of a Diophantine equation.
	/// </summary>
	public struct DiophantineSolution
	{
		/// <summary>
		/// Indicates whether the equation has a solution.
		/// </summary>
		public bool IsSolvable { get; set; }

		/// <summary>
		/// The x-coordinate of the particular solution, or null if no solution exists.
		/// </summary>
		public int? X0 { get; set; }

		/// <summary>
		/// The y-coordinate of the particular solution, or null if no solution exists.
		/// </summary>
		public int? Y0 { get; set; }

		/// <summary>
		/// The coefficient of x in the equation.
		/// </summary>
		public readonly int A;

		/// <summary>
		/// The coefficient of y in the equation.
		/// </summary>
		public readonly int B;

		/// <summary>
		/// The constant term in the equation.
		/// </summary>
		public readonly int C;
		
		/// <summary>D: The greatest common divisor gcd(a,b) = d</summary>
		/// <remarks>Solutions to the Diophantine linear equation exists only if C|D: D is divisible by C.</remarks>
		public readonly int? D;

		/// <summary>
		/// Initializes a new instance of the DiophantineSolution struct.
		/// </summary>
		/// <param name="isSolvable">Indicates whether the equation has a solution.</param>
		/// <param name="x0">The x-coordinate of the particular solution.</param>
		/// <param name="y0">The y-coordinate of the particular solution.</param>
		/// <param name="a">The coefficient of x in the equation.</param>
		/// <param name="b">The coefficient of y in the equation.</param>
		/// <param name="c">The constant term in the equation.</param>
		/// <param name="d">The d= gcd(a,b) greatest common divisor found for a and b</param>
		public DiophantineSolution(bool isSolvable, int? x0, int? y0, int a, int b, int c, int d)
		{
			IsSolvable = isSolvable;
			X0 = x0;
			Y0 = y0;
			A = a;
			B = b;
			C = c;
			D = d;
		}

		/// <summary>
		/// Returns a string representation of the solution.
		/// </summary>
		/// <returns>A string representation of the solution.</returns>
		public override string ToString()
		{
			if (IsSolvable)
			{
				string result = $"""
				Diophantine Linear Equation Solver result for equation 
				AX + BY = C
				{A}x + {B}y = {C} 
				A={A}, B={B}, C={C}
				Greatest common divisor :
				d = gcd(a,b) = {D?.ToString() ?? "No gcd found"}
				Solutions exists? {IsSolvable}, since d|c = {D}|{C}: {C} is divisible by {D}
				Particular solution A(x0) + B(y0) = C
				x0 = {X0}, y0 = {Y0}
				({A}*{X0}) + ({B}*{Y0}) = {C}
				({A * X0}) + ({B * Y0}) = {C}
				Generalized solution:				
				{GetParameterizedSolution()}
				Additional solutions:
				""";
				if (D > 0)
				{
					for (int t = 1; t <= 3; t++)
					{
						int xn = X0!.Value + (B / D!.Value) * t;
						int yn = Y0!.Value - (A / D!.Value) * t;
						result += $"\nSolution {t}: x = {xn}, y = {yn} (t = {t}) , since : ({xn}*{X0}) + ({yn}*{Y0}) = ({xn * X0}) + ({yn * Y0}) = {C}";
					}
				}
				else {
					result += "No additional solutions exists"; 
				}
				return result;
			}
			else
			{
				return $"Diophantine Linear Equation Solver result\nNo solution exists. d = gcd(a,b) = {D}. Solution exists if: {C}|{D}. This gives a remainder of {C % D} != 0. {C} is not divisible by {D}. C = {C} = n*D + r = ({C / D}*{D}) + {C % D}.";
			}
		}

		/// <summary>
		/// Returns the parameterized solution expression as a string.
		/// </summary>
		/// <returns>A string representation of the parameterized solution.</returns>
		public string GetParameterizedSolution()
		{
			if (D == 0){
				return "No solution exists, so there is not general parameterized solution";
			}
			string generalEquation = "x' = x_0 + (B/D), y' = y_0 - (A/D)";
			return $"{generalEquation}{Environment.NewLine}x' = {X0} + ({B}/{D})t, y' = {Y0} - ({A}/{D})t, where t is any integer";
		}
	}

	/// <summary>
	/// Solves the Diophantine equation ax + by = c using the Extended Euclidean Algorithm.
	/// </summary>
	/// <param name="a">The coefficient of x.</param>
	/// <param name="b">The coefficient of y.</param>
	/// <param name="c">The constant term.</param>
	/// <returns>A DiophantineSolution struct containing the solution.</returns>
	public static DiophantineSolution ExtendedGcd(int a, int b, int c)
	{
		int X0, Y0;
		int gcd = ExtendedGcd(a, b, out X0, out Y0);

		if (c % gcd != 0)
		{
			return new DiophantineSolution(false, null, null, a, b, c, gcd);
		}
		else
		{
			int x = X0 * (c / gcd);
			int y = Y0 * (c / gcd);
			return new DiophantineSolution(true, x, y, a, b, c, gcd);
		}
	}

	/// <summary>
	/// Computes the greatest common divisor of a and b, and finds coefficients X0 and Y0 such that a*X0 + b*Y0 = gcd(a, b).
	/// </summary>
	/// <param name="a">The first integer.</param>
	/// <param name="b">The second integer.</param>
	/// <param name="X0">The coefficient of a in the linear combination.</param>
	/// <param="Y0">The coefficient of b in the linear combination.</param>
	/// <returns>The greatest common divisor of a and b.</returns>
	private static int ExtendedGcd(int a, int b, out int X0, out int Y0)
	{
		if (b == 0)
		{
			X0 = 1;
			Y0 = 0;
			return a; // Base case: gcd(a, 0) = a
		}

		int X1, Y1;
		int gcd = ExtendedGcd(b, a % b, out X1, out Y1); // Recursive call

		X0 = Y1; // Back substitution: X0 = Y1
		Y0 = X1 - (a / b) * Y1; // Back substitution: Y0 = X1 - (a / b) * Y1

		return gcd; // Return the gcd
	}
}

The extended Euclidean algorithm calculates the greatest common divisor and in the recursion also applies back substitution. The coefficients X0 and Y0 founds satisfy the equation :

(a · X₀) + (b · Y₀) = gcd(a, b)

That is, a linear combination of the original integers a and b equaling the gcd(a,b). The output of running the code shown above is:

Diophantine Linear Equation Solver result for equation
AX + BY = C
7x + -9y = 3
A=7, B=-9, C=3
Greatest common divisor :
d = gcd(a,b) = 1
Solutions exists? True, since d|c = 1|3: 3 is divisible by 1
Particular solution A(x0) + B(y0) = C
x0 = 12, y0 = 9
(7*12) + (-9*9) = 3
(84) + (-81) = 3
Generalized solution:               
x' = x_0 + (B/D), y' = y_0 - (A/D)
x' = 12 + (-9/1)t, y' = 9 - (7/1)t, where t is any integer
Additional solutions:
Solution 1: x = 3, y = 2 (t = 1) , since : (3*12) + (2*9) = (36) + (18) = 3
Solution 2: x = -6, y = -5 (t = 2) , since : (-6*12) + (-5*9) = (-72) + (-45) = 3
Solution 3: x = -15, y = -12 (t = 3) , since : (-15*12) + (-12*9) = (-180) + (-108) = 3

So an example of a solution of the equation : 7x -9y = 3 Is: x = 12, y = 9, since (7*12) - (9*9) = 84 - 81 = 3. Note that once you have found a solution to a linear Diophantine equation, you can found infinite other solutions. More formally:

Diophantine Equation (Linear) Definition and Solvability

A linear Diophantine equation is an equation of the form:

(a · x) + (b · y) = c

where a, b, and c are integers, and x and y are unknown integers.

Solvability

The linear Diophantine equation (a · x) + (b · y) = c has integer solutions if and only if the greatest common divisor (gcd) of a and b divides c. In mathematical notation:

d is gcd. d|c => gcd(a, b) | c

If gcd(a, b) divides c (it is common to use d to notify gcd as a short term notation), then there exist integer solutions x and y such that:

(a · x) + (b · y) = c

Otherwise, there are no integer solutions.