I wonder if this will lead us down the same path of adding a special flag for array zeroing vs IEEE compliance as with non-complex SCAL :(

You reminded me that we also need to add flags for cscal and zscal.
For the following code:

#include <stdio.h>
#include "cblas.h"
#include <complex.h>
int main() {
    // 向量长度
    int N = 1;

    // 复数缩放因子 alpha
    float alpha[2] = {0.0, 0.0}; 

    // 复数向量 X
    float X[2] = {NAN, NAN};

    // 打印缩放前的向量
    printf("Before scaling:\n");
    for (int i = 0; i < N; i++) {
        printf("X[%d] = %f + %f\n", i, X[i], X[i + 1]);
    }

    // 缩放向量
    cblas_cscal(N, alpha, X, 1);

    // 打印缩放后的向量
    printf("\nAfter scaling:\n");
    for (int i = 0; i < N; i++) {
        printf("X[%d] = %f + %f\n", i, X[i], X[i + 1]);
    }

    return 0;
}

Test output using MKL 2024.2 version：

Before scaling:
X[0] = nan + nan

After scaling:
X[0] = nan + nan

The same output as MKL when using reference BLAS Version 3.12.0.
Test output using OpenBLAS V0.3.29(TARGET=HASWELL)：

Before scaling:
X[0] = nan + nan

After scaling:
X[0] = 0.000000 + 0.000000

This PR will introduce new issues to s/zscal and needs to be revised. (It seems that other platforms also need modifications to avoid the above issues.)

I submitted a PR #5081 attempting to fix the implementation in C.

Thanks, I'll try to take a stab at the other implementations over the weekend.

Unfortunately this appears to have broken most of the pre-existing special handling of 0*NAN and 0*INF cases (which mostly covered only one of the real and imaginary components being "special") ?

I removed some special-case handling code because I felt its correctness was questionable.
Test the following code using MKL 2024.02:

#include <stdio.h>
#include <mkl_cblas.h>
#include <complex.h>
#include <string.h>

int main() {
    int i;
    int N = 17;
    int incX = 1;
    int incY = 1;
    float alpha[2] = {0.0, 0.0};
    float beta[2] = {2.0, 0.0};
    char  trans = 'N';
    float A[17 * 17 * 4];
    float X[17 * 2];
    float Y[17 * 2];

    memset(A, 0, sizeof(A));
    memset(X, 0, sizeof(X));
    for (i = 0; i < (2 * N - 2); i += 4)
    {
        Y[i]     = NAN;
        Y[i + 1] = 1.0;

        Y[i + 2] = INFINITY;
        Y[i + 3] = 1.0;
    }
    Y[2 * N - 2] = NAN;
    Y[2 * N - 1] = 1.0;
    cblas_cgemv(CblasColMajor, CblasNoTrans, N, N, alpha, A, N, X, incX, beta, Y, incY);
    for (i = 0; i < 2 * N; i += 2) {
            printf("%f, %f, ", Y[i], Y[i+1]);
    }
}

The output is:
nan, nan, inf, -nan, nan, nan, inf, -nan, nan, nan, inf, -nan, nan, nan, inf, -nan, nan, nan, inf, -nan, nan, nan, inf, -nan, nan, nan, inf, -nan, nan, nan, inf, -nan, nan, nan
The same output as MKL when using reference BLAS Version 3.12.0.
However, when using OpenBLAS v0.3.29 (TARGET=GENERIC, using the C kernel), the output result is:
nan, 2.000000, inf, 2.000000, nan, 2.000000, inf, 2.000000, nan, 2.000000, inf, 2.000000, nan, 2.000000, inf, 2.000000, nan, 2.000000, inf, 2.000000, nan, 2.000000, inf, 2.000000, nan, 2.000000, inf, 2.000000, nan, 2.000000, inf, 2.000000, nan, 2.000000
I think for the upper-level interfaces that call c/zscal_k (excluding c/zscal), we should also no longer reduce computations.

Adding special value checks for each number is likely to cause a performance drop and make optimization with assembly much more difficult. Should we consider simplifying it?

(spurious merge conflict was caused by pengxu's #5248 which contained your changes to cscal_lasx.c - sorry for not noticing this at the time I merged that PR)

No problem, Just wanted to check if you'd like me to revert the merged changes in my PR.

I think it is all sorted, thank you.