diff --git a/BREAKING-CHANGES.md b/BREAKING-CHANGES.md index 022dc7d2d..4af73b5e6 100644 --- a/BREAKING-CHANGES.md +++ b/BREAKING-CHANGES.md @@ -305,14 +305,16 @@ bracketing for. The change only ever adds `\left(`/`\right)` groups, which CShar parses, so nothing downstream needs a matching change ([#822](https://github.com/asc-community/AngouriMath/issues/822)). -### `Factor` factors a polynomial in two variables +### `Factor` factors a polynomial in more than one variable -After the content is taken out, what remains may still have polynomial coefficients — and where it -is in two variables it can be factored anyway, by **Kronecker's substitution**. A factor of a -polynomial of degree `d` in `x` has degree at most `d` in `x`, so with `s = d + 1` the map -`x^i y^j → t^(i + s*j)` is injective on every monomial that can appear in the polynomial or in any -of its factors. The one-variable image is factored by the existing factoriser, and each subset of -its irreducible factors names a candidate. +After the content is taken out, what remains may still have polynomial coefficients — and it can be +factored anyway, by **Kronecker's substitution written in mixed radix**. A factor of a polynomial +has degree at most `d_i` in each variable `v_i`, because a factor divides it. So with radices +`d_i + 1` and place values `s_0 = 1`, `s_(i+1) = s_i * (d_i + 1)`, the map sending a monomial to +`t^(Σ e_i · s_i)` writes each exponent as one digit of a numeral, and is therefore injective on +every monomial that can appear in the polynomial or in any of its factors. The one-variable image +is factored by the existing factoriser, and each subset of its irreducible factors names a +candidate. | | 2.3.0 | now | |---|---|---| @@ -322,17 +324,23 @@ its irreducible factors names a candidate. | `Factor("x ^ 4 - y ^ 4", "x")` | `null` | `(x + y) * (x ^ 2 + y ^ 2) * (x - y)` | | `Factor("x ^ 2 * y ^ 2 - 1", "x")` | `null` | `(x * y + 1) * (x * y - 1)` | | `Factor("x ^ 2 - y ^ 2 + 2 * x + 1", "x")` | `null` | `(x + y + 1) * (x - y + 1)` | +| `Factor("x ^ 2 - (y + z) ^ 2", "x")` | `null` | `(x + y + z) * (x - y - z)` | +| `Factor("x ^ 2 + 2 * x * y + y ^ 2 - z ^ 2", "x")` | `null` | `(x + y + z) * (x + y - z)` | +| `Factor("(x + y) * (x + z) * (x + w)", "x")` | `null` | `(x + y) * (w + x) * (x + z)` | | `Factor("x ^ 2 + y ^ 2", "x")` | `null` | `null` — irreducible over ℚ | -| `Factor("x * y + z", "x")` | `null` | `null` — three variables | +| `Factor("x * y + z", "x")` | `null` | `null` — irreducible over ℚ | **It cannot answer wrongly.** The substitution is injective on monomials but not on factorisations, so the image may factor further than the polynomial does and a candidate is a guess. Every one is tested by exact division before it is kept, and the assembled factors are divided back into the input, so the failure mode is a refusal. -**What it refuses.** The image has degree `d + s*e` for degree `e` in the second variable, and the -one-variable factoriser stops at 32 — so this reaches bidegrees like (2, 10), (3, 7) and (5, 4) and -refuses past them. The recombination is over subsets, so the image's factor count is capped too. +**What it refuses.** The image has degree `Π (d_i + 1) - 1`, a **product** and not a sum, and the +one-variable factoriser stops at 32 — so the ceiling closes quickly as variables are added. Two +variables reach bidegrees like (2, 10), (3, 7) and (5, 4); three variables of degree 2 fit (27) and +four do not (81). `Factor("x ^ 12 - y ^ 12", "x")` and +`Factor("(x + y + z + w) * (x - y)", "x")` are both `null` for this reason, though both factor +mathematically. The recombination is over subsets, so the image's factor count is capped too. Lifting that ceiling is Hensel lifting with an evaluation homomorphism, which is a different piece of work. @@ -381,14 +389,14 @@ the other variables — is now taken out first, using the same multivariate mach | `Factor("x ^ 2 * y + x * y", "x")` | `null` | `y * x * (x + 1)` | | `Factor("a * x ^ 2 + a * x", "x")` | `null` | `a * x * (x + 1)` | | `Factor("x ^ 2 * y ^ 2 - y ^ 2", "x")` | `null` | `y ^ 2 * (x + 1) * (x - 1)` | -| `Factor("x ^ 2 - y ^ 2", "x")` | `null` | `null` | | `Factor("x * y + z", "x")` | `null` | `null` | -**Only a refusal becomes an answer.** Nothing that already factorised changes, because the new path -runs only where the old one returned `null`. And it is still a refusal wherever the content is a -constant: `x ^ 2 - y ^ 2` genuinely needs factorisation over ℚ(y), which this is not and does not -claim to be. That remains the open half of -[#746](https://github.com/asc-community/AngouriMath/issues/746) item 43. +**Only a refusal becomes an answer.** Nothing that already factorised changes, because this path +runs only where the old one returned `null`. + +Taking the content out does nothing where the content is a constant, so `x ^ 2 - y ^ 2` is not +answered by this change — it needs factorisation over ℚ(y). That is what Kronecker's substitution +does, in the entry above, and the two paths are tried in that order. The test that pinned the refusal carried a comment saying that handing `x * y + y` back *"would say that `y * (x + 1)` does not exist, which is a wrong answer and not a graceful failure"*. It now diff --git a/Sources/.editorconfig b/Sources/.editorconfig index 10473aee5..bc7424bab 100644 --- a/Sources/.editorconfig +++ b/Sources/.editorconfig @@ -27,7 +27,7 @@ file_header_template=\nCopyright (c) 2019-2026 Angouri.\nAngouriMath is licensed [Tests/UnitTests/Core/Transformations/*.cs] file_header_template=\nCopyright (c) 2019-2026 Angouri.\nAngouriMath is licensed under MIT.\nDetails: https://github.com/asc-community/AngouriMath/blob/master/LICENSE.md.\nWebsite: https://am.angouri.org.\n -[AngouriMath/Functions/Algebra/Polynomials/BivariateFactorization.cs] +[AngouriMath/Functions/Algebra/Polynomials/KroneckerFactorization.cs] file_header_template=\nCopyright (c) 2019-2026 Angouri.\nAngouriMath is licensed under MIT.\nDetails: https://github.com/asc-community/AngouriMath/blob/master/LICENSE.md.\nWebsite: https://am.angouri.org.\n [AngouriMath/Functions/Algebra/Groebner/*.cs] diff --git a/Sources/AngouriMath/Convenience/MathS.Polynomials.cs b/Sources/AngouriMath/Convenience/MathS.Polynomials.cs index daa0af086..4b2c017b6 100644 --- a/Sources/AngouriMath/Convenience/MathS.Polynomials.cs +++ b/Sources/AngouriMath/Convenience/MathS.Polynomials.cs @@ -133,10 +133,10 @@ public static class Polynomials /// So the content in — the greatest common divisor of the /// coefficients, which is a polynomial in the other variables — is taken out first, /// using the same multivariate machinery is built from, and what - /// remains goes down the ordinary path. Where the content is a constant this has - /// nothing to offer and says so, which is the honest answer for - /// x ^ 2 - y ^ 2: that one genuinely needs factorisation over ℚ(y) and is not - /// what this does. + /// remains goes down the ordinary path. Where the content is a constant that path + /// has nothing to offer, and answers instead — + /// x ^ 2 - y ^ 2 is (x + y) * (x - y), which is a factorisation over + /// ℚ(y) reached by substitution rather than by lifting. /// /// private static Entity? FactorAfterTakingOutTheContent(Entity expr, Variable variable) @@ -157,12 +157,12 @@ public static class Polynomials if (poly.DivideExact(content) is not { } primitive) return null; - // What is left may still have polynomial coefficients, and where it is in two - // variables it can be factored anyway -- see BivariateFactorization. + // What is left may still have polynomial coefficients, and it can be factored + // anyway while the substitution's ceiling allows -- see KroneckerFactorization. var rest = Assemble( PolynomialFactorization.FactorComplete(primitive.ToEntity(variables), variable), variable) - ?? Bivariate(primitive, variables, index, variable); + ?? Kronecker(primitive, variables, index, variable); if (rest is null) return null; return content.IsConstant && content.DivideExact(content) is not null @@ -177,21 +177,19 @@ private static bool SameAsOne(MultivariatePolynomial poly) => poly.IsConstant && poly.CoefficientOf(0).CompareTo(ERational.One) == 0; /// - /// The factorisation of a polynomial in exactly two variables, as an expression. + /// The factorisation of a polynomial in more than one variable, as an expression. /// /// - /// Kronecker's substitution: see for what it + /// Kronecker's substitution: see for what it /// does, what it refuses, and why a wrong answer is not among the things it can do. /// - private static Entity? Bivariate( + private static Entity? Kronecker( MultivariatePolynomial poly, IReadOnlyList variables, IReadOnlyDictionary index, Variable variable) { - if (variables.Count != 2) + if (variables.Count < 2) return null; - var main = index[variable]; - var other = main == 0 ? 1 : 0; - if (BivariateFactorization.Factor(poly, main, other) is not { } factors + if (KroneckerFactorization.Factor(poly, index[variable]) is not { } factors || factors.Count < 2) return null; // Repeated factors are collected into a power, as the one-variable path does: diff --git a/Sources/AngouriMath/Functions/Algebra/Polynomials/BivariateFactorization.cs b/Sources/AngouriMath/Functions/Algebra/Polynomials/KroneckerFactorization.cs similarity index 55% rename from Sources/AngouriMath/Functions/Algebra/Polynomials/BivariateFactorization.cs rename to Sources/AngouriMath/Functions/Algebra/Polynomials/KroneckerFactorization.cs index dc8b8bd48..9066bfc31 100644 --- a/Sources/AngouriMath/Functions/Algebra/Polynomials/BivariateFactorization.cs +++ b/Sources/AngouriMath/Functions/Algebra/Polynomials/KroneckerFactorization.cs @@ -11,38 +11,41 @@ namespace AngouriMath.Functions { /// - /// Factorisation of a polynomial in two variables, by reducing it to one variable and - /// putting the answer back. + /// Factorisation of a polynomial in any number of variables, by reducing it to one variable + /// and putting the answer back. /// /// /// - /// Kronecker's substitution. A factor of a polynomial of degree d in x - /// has degree at most d in x, so with s = d + 1 the map - /// x^i y^j -> t^(i + s*j) is injective on every monomial that can appear in the - /// polynomial or in any of its factors: i is the remainder and j the quotient of - /// the exponent by s, and neither can be confused with another pair. So a factorisation - /// of the one-variable image can be read back, and each subset of its irreducible factors - /// names a candidate. + /// Kronecker's substitution, in mixed radix. A factor of a polynomial has degree at + /// most d_i in each variable v_i, because a factor divides it. So with radices + /// d_i + 1 and place values s_0 = 1, s_(i+1) = s_i * (d_i + 1), the map + /// v_0^e_0 · … · v_(k-1)^e_(k-1) -> t^(Σ e_i · s_i) writes each exponent as one digit + /// of a numeral and is therefore injective on every monomial that can appear in the + /// polynomial or in any of its factors. A factorisation of the one-variable image can be read + /// back digit by digit, and each subset of its irreducible factors names a candidate. /// /// - /// A candidate is a guess and is checked by division. The image can factor further than - /// the polynomial does — the substitution is injective on monomials, not on factorisations — - /// so a subset whose product reads back as a polynomial need not divide the original. Every - /// one is tested with exact division before it is kept, which is why this cannot answer - /// wrongly: the worst it does is fail to find a factorisation that exists and say so. + /// A candidate is a guess and is checked by division. The image can factor further + /// than the polynomial does — the substitution is injective on monomials, not on + /// factorisations — so a subset whose product reads back as a polynomial need not divide the + /// original. Every one is tested with exact division before it is kept, which is why this + /// cannot answer wrongly: the worst it does is fail to find a factorisation that exists and + /// say so. /// /// - /// What it will not do. The image has degree d + s * e for a polynomial of - /// degree e in y, and the one-variable factoriser stops at - /// — so this reaches bidegrees like (2, 10), - /// (3, 7) and (5, 4) and refuses beyond them. The recombination is over subsets, so the - /// number of irreducible factors of the image is capped as well. Both limits are refusals, - /// never wrong answers, and neither is the algorithm one would write to lift this ceiling: - /// that is Hensel lifting with an evaluation homomorphism, and it is a different piece of - /// work (#746 item 43). + /// What it will not do. The image has degree Π (d_i + 1) - 1, a *product* and + /// not a sum, and the one-variable factoriser stops at + /// — so the ceiling closes quickly as variables are + /// added. Two variables reach bidegrees like (2, 10), (3, 7) and (5, 4); three variables of + /// degree 2 fit (27 ≤ 32) and four do not (81); and a quadratic in eight variables is far + /// past it. The recombination is over subsets, so the number of irreducible factors of the + /// image is capped as well. Both limits are refusals, never wrong answers, and neither is the + /// algorithm one would write to lift them: that is Hensel lifting with an evaluation + /// homomorphism, and it is a different piece of work + /// (#746 item 43). /// /// - internal static class BivariateFactorization + internal static class KroneckerFactorization { /// /// Beyond this the subset search is refused rather than paid for: the recombination is @@ -52,31 +55,47 @@ internal static class BivariateFactorization private const int MaxImageFactors = 12; /// - /// The factors of in and , - /// each to the first power and each of positive degree in , or - /// where nothing could be settled. A single factor means the - /// polynomial did not factor, which is an answer. + /// The factors of in and whatever other + /// variables it has, each to the first power and each of positive degree in + /// , or where nothing could be settled. A + /// single factor means the polynomial did not factor, which is an answer. /// internal static IReadOnlyList? Factor( - MultivariatePolynomial poly, int x, int y) + MultivariatePolynomial poly, int main) { - var degreeInX = poly.DegreeIn(x); - var degreeInY = poly.DegreeIn(y); - if (poly.IsZero || degreeInX < 1 || degreeInY < 1) + if (poly.IsZero || poly.DegreeIn(main) < 1) return null; - var stride = degreeInX + 1; - var imageDegree = degreeInX + stride * degreeInY; - if (imageDegree > IntegerPolynomial.MaxDegree) + // The main variable is placed first so that its exponent is the lowest digit; the + // rest follow in index order, and a variable the polynomial does not use is left out + // rather than given a radix of one. + var order = new List { main }; + for (var variable = 0; variable < poly.VariableCount; variable++) + if (variable != main && poly.DegreeIn(variable) > 0) + order.Add(variable); + if (order.Count < 2) return null; - if (ToImage(poly, x, y, stride, imageDegree) is not { } image) + var radices = new int[order.Count]; + var places = new int[order.Count]; + long size = 1; + for (var i = 0; i < order.Count; i++) + { + radices[i] = poly.DegreeIn(order[i]) + 1; + places[i] = (int)size; + size *= radices[i]; + if (size > IntegerPolynomial.MaxDegree + 1) + return null; + } + var imageDegree = (int)size - 1; + + if (ToImage(poly, order, places, radices, imageDegree) is not { } image) return null; if (PolynomialFactorization.FactorPrimitive(image.PrimitivePart()) is not { } parts) return null; // The multiplicities are flattened: a square in the image may or may not be a square - // in two variables, and the recombination settles that by division rather than by + // in the original, and the recombination settles that by division rather than by // carrying the exponent across the substitution. var irreducibles = new List(); foreach (var part in parts) @@ -89,7 +108,7 @@ internal static class BivariateFactorization if (irreducibles.Count < 2) return new[] { poly }; - return Recombine(poly, irreducibles, x, y, stride); + return Recombine(poly, irreducibles, main, order, places, radices); } /// @@ -97,23 +116,15 @@ internal static class BivariateFactorization /// not wanted, since a rational multiple of a factor is the same factor. /// private static IntegerPolynomial? ToImage( - MultivariatePolynomial poly, int x, int y, int stride, int imageDegree) + MultivariatePolynomial poly, IReadOnlyList order, + IReadOnlyList places, IReadOnlyList radices, int imageDegree) { var coefficients = new ERational[imageDegree + 1]; for (var i = 0; i < coefficients.Length; i++) coefficients[i] = ERational.Zero; - foreach (var byX in poly.CoefficientsIn(x)) - foreach (var byY in byX.Value.CoefficientsIn(y)) - { - // Anything left is a third variable, and this is the two-variable case. - if (!byY.Value.IsConstant) - return null; - var at = byX.Key + stride * byY.Key; - if (at > imageDegree) - return null; - coefficients[at] = coefficients[at].Add(byY.Value.CoefficientOf(0)); - } + if (!Peel(poly, 0, 0)) + return null; var denominator = EInteger.One; foreach (var coefficient in coefficients) @@ -123,32 +134,57 @@ internal static class BivariateFactorization whole[i] = coefficients[i].Numerator .Multiply(denominator.Divide(coefficients[i].Denominator)); return IntegerPolynomial.Create(whole); + + // One variable at a time, adding that exponent's digit to the place already + // accumulated. What is left when every variable has been peeled has to be a constant. + bool Peel(MultivariatePolynomial rest, int depth, int at) + { + if (depth == order.Count) + { + if (!rest.IsConstant) + return false; + coefficients[at] = coefficients[at].Add(rest.CoefficientOf(0)); + return true; + } + foreach (var pair in rest.CoefficientsIn(order[depth])) + { + if (pair.Key >= radices[depth]) + return false; + if (!Peel(pair.Value, depth + 1, at + places[depth] * pair.Key)) + return false; + } + return true; + } } private static EInteger Lcm(EInteger left, EInteger right) => left.Divide(left.Gcd(right)).Multiply(right); /// - /// Two variables again, reading each exponent as its remainder and quotient by - /// — or where a power is past what a - /// packed monomial holds. + /// Reading each exponent back as one digit of the numeral — or + /// where a power is past what a packed monomial holds. /// private static MultivariatePolynomial? FromImage( - IntegerPolynomial image, int variableCount, int x, int y, int stride) + IntegerPolynomial image, int variableCount, IReadOnlyList order, + IReadOnlyList places, IReadOnlyList radices) { var result = MultivariatePolynomial.Zero(variableCount); for (var power = 0; power <= image.Degree; power++) { if (image[power].IsZero) continue; - var inX = power % stride; - var inY = power / stride; - if (inX > MultivariatePolynomial.MaxDegree || inY > MultivariatePolynomial.MaxDegree) - return null; - if (MultivariatePolynomial.Monomial(variableCount, x).Power(inX) is not { } partX - || MultivariatePolynomial.Monomial(variableCount, y).Power(inY) is not { } partY - || partX.Multiply(partY) is not { } monomial) - return null; + var monomial = MultivariatePolynomial.One(variableCount); + for (var i = 0; i < order.Count; i++) + { + var digit = power / places[i] % radices[i]; + if (digit > MultivariatePolynomial.MaxDegree) + return null; + if (MultivariatePolynomial.Monomial(variableCount, order[i]).Power(digit) + is not { } part + || monomial.Multiply(part) is not { } multiplied) + return null; + monomial = multiplied; + } result = result.Add(monomial.ScaleBy(ERational.Create(image[power], EInteger.One))); } return result.IsZero ? null : result; @@ -164,7 +200,8 @@ private static EInteger Lcm(EInteger left, EInteger right) /// after each success because the remaining polynomial has changed. /// private static IReadOnlyList? Recombine( - MultivariatePolynomial poly, List irreducibles, int x, int y, int stride) + MultivariatePolynomial poly, List irreducibles, int main, + IReadOnlyList order, IReadOnlyList places, IReadOnlyList radices) { var found = new List(); var remaining = poly; @@ -183,8 +220,9 @@ private static EInteger Lcm(EInteger left, EInteger right) product = multiplied; else return null; - if (FromImage(product, poly.VariableCount, x, y, stride) is not { } candidate - || candidate.DegreeIn(x) < 1 + if (FromImage(product, poly.VariableCount, order, places, radices) + is not { } candidate + || candidate.DegreeIn(main) < 1 || remaining.DivideExact(candidate) is not { } quotient) continue; found.Add(candidate.Normalized()); diff --git a/Sources/Tests/UnitTests/Algebra/Polynomials/PolynomialSurfaceTest.cs b/Sources/Tests/UnitTests/Algebra/Polynomials/PolynomialSurfaceTest.cs index ac7d2f853..ae0bb895e 100644 --- a/Sources/Tests/UnitTests/Algebra/Polynomials/PolynomialSurfaceTest.cs +++ b/Sources/Tests/UnitTests/Algebra/Polynomials/PolynomialSurfaceTest.cs @@ -157,15 +157,19 @@ public void ASquareFreePartHasNoRepeatedRoot(string input) } /// - /// A polynomial in two variables is factored by Kronecker's substitution: with - /// s one more than its degree in x, the map x^i y^j -> t^(i + s*j) - /// is injective on every monomial that can appear in it or in any of its factors, so a - /// factorisation of the image reads back and each subset of its irreducible factors - /// names a candidate. Every candidate is checked by exact division, so the failure mode - /// is a refusal and not a wrong answer. + /// A polynomial in several variables is factored by Kronecker's substitution, written in + /// mixed radix: with radices d_i + 1 and place values s_0 = 1, + /// s_(i+1) = s_i * (d_i + 1), the map sending a monomial to + /// t^(sum of e_i * s_i) writes each exponent as one digit of a numeral, so it is + /// injective on every monomial that can appear in the polynomial or in any of its + /// factors. A factorisation of the one-variable image reads back digit by digit, and each + /// subset of its irreducible factors names a candidate. Every candidate is checked by + /// exact division, so the failure mode is a refusal and not a wrong answer. /// /// - /// x ^ 2 - y ^ 2 is the case #746 item 43 names, and the one this test exists for. + /// x ^ 2 - y ^ 2 is the case #746 item 43 names. The three- and four-variable rows + /// are the generalisation: the exponent vector is a numeral whatever its length, and only + /// the image's degree — a product of the radices, not a sum — decides what fits. /// Compared numerically, for the reason the content test above gives. /// [Theory] @@ -175,7 +179,12 @@ public void ASquareFreePartHasNoRepeatedRoot(string input) [InlineData("x ^ 2 * y ^ 2 - 1", 2)] [InlineData("x ^ 4 - y ^ 4", 3)] [InlineData("x ^ 2 - y ^ 2 + 2 * x + 1", 2)] - public void APolynomialInTwoVariablesIsFactored(string input, int distinctFactors) + [InlineData("x ^ 2 - (y + z) ^ 2", 2)] + [InlineData("x ^ 2 + 2 * x * y + y ^ 2 - z ^ 2", 2)] + [InlineData("x * y - x - y + 1", 2)] + [InlineData("(x + y) * (x + z + w)", 2)] + [InlineData("(x + y) * (x + z) * (x + w)", 3)] + public void APolynomialInSeveralVariablesIsFactored(string input, int distinctFactors) { var expr = input.ToEntity(); var factored = MathS.Polynomials.Factor(expr, "x"); @@ -207,6 +216,38 @@ public void APolynomialInTwoVariablesIsFactored(string input, int distinctFactor private static int CountFactors(Entity product) => product is Mulf(var left, var right) ? CountFactors(left) + CountFactors(right) : 1; + /// + /// The image's degree is a product of the radices and not a sum, so the ceiling + /// closes quickly as variables are added — and where it closes the answer is a refusal. + /// + /// + /// Every row here does factor mathematically and is declined anyway, which is the shape + /// of every limit in this path: a refusal is a possible answer and a wrong one is not. + /// The degrees are (2, 2, 1, 1) and (4, 1, 1, 1, 1), giving images of + /// degree 35 and 79 against an IntegerPolynomial.MaxDegree of 32; the third is + /// two variables and past it on its own. + /// + [Theory] + [InlineData("(x + y + z + w) * (x - y)")] + [InlineData("(x + y) * (x + z) * (x + w) * (x + v)")] + [InlineData("x ^ 12 - y ^ 12")] + public void PastTheSubstitutionsCeilingItRefuses(string input) + => Assert.Null(MathS.Polynomials.Factor(input.ToEntity(), "x")); + + /// + /// The main variable is a parameter and not a convention: the same polynomial factored + /// with respect to another variable is the same factorisation. + /// + [Fact] + public void TheMainVariableIsAParameter() + { + var inX = MathS.Polynomials.Factor("x ^ 2 - y ^ 2".ToEntity(), "x"); + var inY = MathS.Polynomials.Factor("x ^ 2 - y ^ 2".ToEntity(), "y"); + Assert.NotNull(inX); + Assert.NotNull(inY); + Assert.Equal(inX!.Expand().Simplify(), inY!.Expand().Simplify()); + } + /// /// The square-free part is p / gcd(p, dp/dx) whatever ring the coefficients live /// in, so it is not univariate for any reason but the representation it used to be