Skip to content

Commit ea2eb00

Browse files
tbaederrtbaederr
committed
[clang][Interp] Implement Complex division
1 parent b89efd9 commit ea2eb00

File tree

6 files changed

+165
-33
lines changed

6 files changed

+165
-33
lines changed

‎clang/lib/AST/ExprConstShared.h

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -62,5 +62,8 @@ GCCTypeClass EvaluateBuiltinClassifyType(QualType T,
6262
void HandleComplexComplexMul(llvm::APFloat A, llvm::APFloat B, llvm::APFloat C,
6363
llvm::APFloat D, llvm::APFloat &ResR,
6464
llvm::APFloat &ResI);
65+
void HandleComplexComplexDiv(llvm::APFloat A, llvm::APFloat B, llvm::APFloat C,
66+
llvm::APFloat D, llvm::APFloat &ResR,
67+
llvm::APFloat &ResI);
6568

6669
#endif

‎clang/lib/AST/ExprConstant.cpp

Lines changed: 43 additions & 33 deletions
Original file line numberDiff line numberDiff line change
@@ -15182,6 +15182,48 @@ void HandleComplexComplexMul(APFloat A, APFloat B, APFloat C, APFloat D,
1518215182
}
1518315183
}
1518415184

15185+
void HandleComplexComplexDiv(APFloat A, APFloat B, APFloat C, APFloat D,
15186+
APFloat &ResR, APFloat &ResI) {
15187+
// This is an implementation of complex division according to the
15188+
// constraints laid out in C11 Annex G. The implementation uses the
15189+
// following naming scheme:
15190+
// (a + ib) / (c + id)
15191+
15192+
int DenomLogB = 0;
15193+
APFloat MaxCD = maxnum(abs(C), abs(D));
15194+
if (MaxCD.isFinite()) {
15195+
DenomLogB = ilogb(MaxCD);
15196+
C = scalbn(C, -DenomLogB, APFloat::rmNearestTiesToEven);
15197+
D = scalbn(D, -DenomLogB, APFloat::rmNearestTiesToEven);
15198+
}
15199+
APFloat Denom = C * C + D * D;
15200+
ResR =
15201+
scalbn((A * C + B * D) / Denom, -DenomLogB, APFloat::rmNearestTiesToEven);
15202+
ResI =
15203+
scalbn((B * C - A * D) / Denom, -DenomLogB, APFloat::rmNearestTiesToEven);
15204+
if (ResR.isNaN() && ResI.isNaN()) {
15205+
if (Denom.isPosZero() && (!A.isNaN() || !B.isNaN())) {
15206+
ResR = APFloat::getInf(ResR.getSemantics(), C.isNegative()) * A;
15207+
ResI = APFloat::getInf(ResR.getSemantics(), C.isNegative()) * B;
15208+
} else if ((A.isInfinity() || B.isInfinity()) && C.isFinite() &&
15209+
D.isFinite()) {
15210+
A = APFloat::copySign(APFloat(A.getSemantics(), A.isInfinity() ? 1 : 0),
15211+
A);
15212+
B = APFloat::copySign(APFloat(B.getSemantics(), B.isInfinity() ? 1 : 0),
15213+
B);
15214+
ResR = APFloat::getInf(ResR.getSemantics()) * (A * C + B * D);
15215+
ResI = APFloat::getInf(ResI.getSemantics()) * (B * C - A * D);
15216+
} else if (MaxCD.isInfinity() && A.isFinite() && B.isFinite()) {
15217+
C = APFloat::copySign(APFloat(C.getSemantics(), C.isInfinity() ? 1 : 0),
15218+
C);
15219+
D = APFloat::copySign(APFloat(D.getSemantics(), D.isInfinity() ? 1 : 0),
15220+
D);
15221+
ResR = APFloat::getZero(ResR.getSemantics()) * (A * C + B * D);
15222+
ResI = APFloat::getZero(ResI.getSemantics()) * (B * C - A * D);
15223+
}
15224+
}
15225+
}
15226+
1518515227
bool ComplexExprEvaluator::VisitBinaryOperator(const BinaryOperator *E) {
1518615228
if (E->isPtrMemOp() || E->isAssignmentOp() || E->getOpcode() == BO_Comma)
1518715229
return ExprEvaluatorBaseTy::VisitBinaryOperator(E);
@@ -15319,39 +15361,7 @@ bool ComplexExprEvaluator::VisitBinaryOperator(const BinaryOperator *E) {
1531915361
// No real optimizations we can do here, stub out with zero.
1532015362
B = APFloat::getZero(A.getSemantics());
1532115363
}
15322-
int DenomLogB = 0;
15323-
APFloat MaxCD = maxnum(abs(C), abs(D));
15324-
if (MaxCD.isFinite()) {
15325-
DenomLogB = ilogb(MaxCD);
15326-
C = scalbn(C, -DenomLogB, APFloat::rmNearestTiesToEven);
15327-
D = scalbn(D, -DenomLogB, APFloat::rmNearestTiesToEven);
15328-
}
15329-
APFloat Denom = C * C + D * D;
15330-
ResR = scalbn((A * C + B * D) / Denom, -DenomLogB,
15331-
APFloat::rmNearestTiesToEven);
15332-
ResI = scalbn((B * C - A * D) / Denom, -DenomLogB,
15333-
APFloat::rmNearestTiesToEven);
15334-
if (ResR.isNaN() && ResI.isNaN()) {
15335-
if (Denom.isPosZero() && (!A.isNaN() || !B.isNaN())) {
15336-
ResR = APFloat::getInf(ResR.getSemantics(), C.isNegative()) * A;
15337-
ResI = APFloat::getInf(ResR.getSemantics(), C.isNegative()) * B;
15338-
} else if ((A.isInfinity() || B.isInfinity()) && C.isFinite() &&
15339-
D.isFinite()) {
15340-
A = APFloat::copySign(
15341-
APFloat(A.getSemantics(), A.isInfinity() ? 1 : 0), A);
15342-
B = APFloat::copySign(
15343-
APFloat(B.getSemantics(), B.isInfinity() ? 1 : 0), B);
15344-
ResR = APFloat::getInf(ResR.getSemantics()) * (A * C + B * D);
15345-
ResI = APFloat::getInf(ResI.getSemantics()) * (B * C - A * D);
15346-
} else if (MaxCD.isInfinity() && A.isFinite() && B.isFinite()) {
15347-
C = APFloat::copySign(
15348-
APFloat(C.getSemantics(), C.isInfinity() ? 1 : 0), C);
15349-
D = APFloat::copySign(
15350-
APFloat(D.getSemantics(), D.isInfinity() ? 1 : 0), D);
15351-
ResR = APFloat::getZero(ResR.getSemantics()) * (A * C + B * D);
15352-
ResI = APFloat::getZero(ResI.getSemantics()) * (B * C - A * D);
15353-
}
15354-
}
15364+
HandleComplexComplexDiv(A, B, C, D, ResR, ResI);
1535515365
}
1535615366
} else {
1535715367
if (RHS.getComplexIntReal() == 0 && RHS.getComplexIntImag() == 0)

‎clang/lib/AST/Interp/ByteCodeExprGen.cpp

Lines changed: 28 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -870,6 +870,19 @@ bool ByteCodeExprGen<Emitter>::VisitComplexBinOp(const BinaryOperator *E) {
870870
return this->emitMulc(ElemT, E);
871871
}
872872

873+
if (Op == BO_Div && LHSType->isAnyComplexType() &&
874+
RHSType->isAnyComplexType()) {
875+
assert(classifyPrim(LHSType->getAs<ComplexType>()->getElementType()) ==
876+
classifyPrim(RHSType->getAs<ComplexType>()->getElementType()));
877+
PrimType ElemT =
878+
classifyPrim(LHSType->getAs<ComplexType>()->getElementType());
879+
if (!this->visit(LHS))
880+
return false;
881+
if (!this->visit(RHS))
882+
return false;
883+
return this->emitDivc(ElemT, E);
884+
}
885+
873886
// Evaluate LHS and save value to LHSOffset.
874887
bool LHSIsComplex;
875888
unsigned LHSOffset;
@@ -981,6 +994,21 @@ bool ByteCodeExprGen<Emitter>::VisitComplexBinOp(const BinaryOperator *E) {
981994
return false;
982995
}
983996
break;
997+
case BO_Div:
998+
if (!loadComplexValue(LHSIsComplex, false, ElemIndex, LHSOffset, LHS))
999+
return false;
1000+
1001+
if (!loadComplexValue(RHSIsComplex, false, ElemIndex, RHSOffset, RHS))
1002+
return false;
1003+
1004+
if (ResultElemT == PT_Float) {
1005+
if (!this->emitDivf(getRoundingMode(E), E))
1006+
return false;
1007+
} else {
1008+
if (!this->emitDiv(ResultElemT, E))
1009+
return false;
1010+
}
1011+
break;
9841012

9851013
default:
9861014
return false;

‎clang/lib/AST/Interp/Interp.h

Lines changed: 72 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -425,6 +425,78 @@ inline bool Mulc(InterpState &S, CodePtr OpPC) {
425425
return true;
426426
}
427427

428+
template <PrimType Name, class T = typename PrimConv<Name>::T>
429+
inline bool Divc(InterpState &S, CodePtr OpPC) {
430+
const Pointer &RHS = S.Stk.pop<Pointer>();
431+
const Pointer &LHS = S.Stk.pop<Pointer>();
432+
const Pointer &Result = S.Stk.peek<Pointer>();
433+
434+
if constexpr (std::is_same_v<T, Floating>) {
435+
APFloat A = LHS.atIndex(0).deref<Floating>().getAPFloat();
436+
APFloat B = LHS.atIndex(1).deref<Floating>().getAPFloat();
437+
APFloat C = RHS.atIndex(0).deref<Floating>().getAPFloat();
438+
APFloat D = RHS.atIndex(1).deref<Floating>().getAPFloat();
439+
440+
APFloat ResR(A.getSemantics());
441+
APFloat ResI(A.getSemantics());
442+
HandleComplexComplexDiv(A, B, C, D, ResR, ResI);
443+
444+
// Copy into the result.
445+
Result.atIndex(0).deref<Floating>() = Floating(ResR);
446+
Result.atIndex(0).initialize();
447+
Result.atIndex(1).deref<Floating>() = Floating(ResI);
448+
Result.atIndex(1).initialize();
449+
Result.initialize();
450+
} else {
451+
// Integer element type.
452+
const T &LHSR = LHS.atIndex(0).deref<T>();
453+
const T &LHSI = LHS.atIndex(1).deref<T>();
454+
const T &RHSR = RHS.atIndex(0).deref<T>();
455+
const T &RHSI = RHS.atIndex(1).deref<T>();
456+
unsigned Bits = LHSR.bitWidth();
457+
const T Zero = T::from(0, Bits);
458+
459+
if (Compare(RHSR, Zero) == ComparisonCategoryResult::Equal &&
460+
Compare(RHSI, Zero) == ComparisonCategoryResult::Equal) {
461+
const SourceInfo &E = S.Current->getSource(OpPC);
462+
S.FFDiag(E, diag::note_expr_divide_by_zero);
463+
return false;
464+
}
465+
466+
// Den = real(RHS)² + imag(RHS)²
467+
T A, B;
468+
if (T::mul(RHSR, RHSR, Bits, &A) || T::mul(RHSI, RHSI, Bits, &B))
469+
return false;
470+
T Den;
471+
if (T::add(A, B, Bits, &Den))
472+
return false;
473+
474+
// real(Result) = ((real(LHS) * real(RHS)) + (imag(LHS) * imag(RHS))) / Den
475+
T &ResultR = Result.atIndex(0).deref<T>();
476+
T &ResultI = Result.atIndex(1).deref<T>();
477+
478+
if (T::mul(LHSR, RHSR, Bits, &A) || T::mul(LHSI, RHSI, Bits, &B))
479+
return false;
480+
if (T::add(A, B, Bits, &ResultR))
481+
return false;
482+
if (T::div(ResultR, Den, Bits, &ResultR))
483+
return false;
484+
Result.atIndex(0).initialize();
485+
486+
// imag(Result) = ((imag(LHS) * real(RHS)) - (real(LHS) * imag(RHS))) / Den
487+
if (T::mul(LHSI, RHSR, Bits, &A) || T::mul(LHSR, RHSI, Bits, &B))
488+
return false;
489+
if (T::sub(A, B, Bits, &ResultI))
490+
return false;
491+
if (T::div(ResultI, Den, Bits, &ResultI))
492+
return false;
493+
Result.atIndex(1).initialize();
494+
Result.initialize();
495+
}
496+
497+
return true;
498+
}
499+
428500
/// 1) Pops the RHS from the stack.
429501
/// 2) Pops the LHS from the stack.
430502
/// 3) Pushes 'LHS & RHS' on the stack

‎clang/lib/AST/Interp/Opcodes.td

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -544,6 +544,10 @@ def Mulc : Opcode {
544544
def Rem : IntegerOpcode;
545545
def Div : IntegerOpcode;
546546
def Divf : FloatOpcode;
547+
def Divc : Opcode {
548+
let Types = [NumberTypeClass];
549+
let HasGroup = 1;
550+
}
547551

548552
def BitAnd : IntegerOpcode;
549553
def BitOr : IntegerOpcode;

‎clang/test/AST/Interp/complex.cpp

Lines changed: 15 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -40,6 +40,21 @@ constexpr _Complex int IIMC = IIMA * IIMB;
4040
static_assert(__real(IIMC) == -30, "");
4141
static_assert(__imag(IIMC) == 40, "");
4242

43+
static_assert(1.0j / 0.0 == 1); // both-error {{static assertion}} \
44+
// both-note {{division by zero}}
45+
static_assert(__builtin_isinf_sign(__real__((1.0 + 1.0j) / (0.0 + 0.0j))) == 1);
46+
static_assert(__builtin_isinf_sign(__real__((1.0 + 1.0j) / 0.0)) == 1); // both-error {{static assertion}} \
47+
// both-note {{division by zero}}
48+
static_assert(__builtin_isinf_sign(__real__((__builtin_inf() + 1.0j) / (0.0 + 0.0j))) == 1);
49+
static_assert(__builtin_isinf_sign(__imag__((1.0 + InfC) / (0.0 + 0.0j))) == 1);
50+
static_assert(__builtin_isinf_sign(__imag__((InfInf) / (0.0 + 0.0j))) == 1);
51+
52+
constexpr _Complex int IIDA = {10,20};
53+
constexpr _Complex int IIDB = {1,2};
54+
constexpr _Complex int IIDC = IIDA / IIDB;
55+
static_assert(__real(IIDC) == 10, "");
56+
static_assert(__imag(IIDC) == 0, "");
57+
4358
constexpr _Complex int Comma1 = {1, 2};
4459
constexpr _Complex int Comma2 = (0, Comma1);
4560
static_assert(Comma1 == Comma1, "");

0 commit comments

Comments
 (0)