diff options
Diffstat (limited to 'packages/base/src/C/lapack-aux.c')
-rw-r--r-- | packages/base/src/C/lapack-aux.c | 177 |
1 files changed, 177 insertions, 0 deletions
diff --git a/packages/base/src/C/lapack-aux.c b/packages/base/src/C/lapack-aux.c index e5e45ef..d56d466 100644 --- a/packages/base/src/C/lapack-aux.c +++ b/packages/base/src/C/lapack-aux.c | |||
@@ -1287,6 +1287,29 @@ int multiplyQ(int ta, int tb, KQMAT(a),KQMAT(b),QMAT(r)) { | |||
1287 | OK | 1287 | OK |
1288 | } | 1288 | } |
1289 | 1289 | ||
1290 | int multiplyI(int ta, int tb, KIMAT(a), KIMAT(b), IMAT(r)) { | ||
1291 | int i,j,k; | ||
1292 | int n; | ||
1293 | int u, v; | ||
1294 | if (ta==0) { | ||
1295 | n = ac; | ||
1296 | } else { | ||
1297 | n = ar; | ||
1298 | } | ||
1299 | for (i=0;i<rr;i++) { | ||
1300 | for (j=0;j<rc;j++) { | ||
1301 | rp[i*rc+j] = 0; | ||
1302 | for (k=0; k<n; k++) { | ||
1303 | u = ta==0 ? ap[i*ac+k] : ap[k*ac+i]; | ||
1304 | v = tb==0 ? bp[k*bc+j] : bp[j*bc+k]; | ||
1305 | rp[i*rc+j] += u*v; | ||
1306 | } | ||
1307 | } | ||
1308 | } | ||
1309 | OK | ||
1310 | } | ||
1311 | |||
1312 | |||
1290 | //////////////////// transpose ///////////////////////// | 1313 | //////////////////// transpose ///////////////////////// |
1291 | 1314 | ||
1292 | int transF(KFMAT(x),FMAT(t)) { | 1315 | int transF(KFMAT(x),FMAT(t)) { |
@@ -1350,6 +1373,19 @@ int transP(KPMAT(x), PMAT(t)) { | |||
1350 | OK | 1373 | OK |
1351 | } | 1374 | } |
1352 | 1375 | ||
1376 | int transI(KIMAT(x),IMAT(t)) { | ||
1377 | REQUIRES(xr==tc && xc==tr,BAD_SIZE); | ||
1378 | DEBUGMSG("transI"); | ||
1379 | int i,j; | ||
1380 | for (i=0; i<tr; i++) { | ||
1381 | for (j=0; j<tc; j++) { | ||
1382 | tp[i*tc+j] = xp[j*xc+i]; | ||
1383 | } | ||
1384 | } | ||
1385 | OK | ||
1386 | } | ||
1387 | |||
1388 | |||
1353 | //////////////////// constant ///////////////////////// | 1389 | //////////////////// constant ///////////////////////// |
1354 | 1390 | ||
1355 | int constantF(float * pval, FVEC(r)) { | 1391 | int constantF(float * pval, FVEC(r)) { |
@@ -1401,6 +1437,18 @@ int constantP(void* pval, PVEC(r)) { | |||
1401 | OK | 1437 | OK |
1402 | } | 1438 | } |
1403 | 1439 | ||
1440 | |||
1441 | int constantI(int * pval, IVEC(r)) { | ||
1442 | DEBUGMSG("constantI") | ||
1443 | int k; | ||
1444 | int val = *pval; | ||
1445 | for(k=0;k<rn;k++) { | ||
1446 | rp[k]=val; | ||
1447 | } | ||
1448 | OK | ||
1449 | } | ||
1450 | |||
1451 | |||
1404 | //////////////////// float-double conversion ///////////////////////// | 1452 | //////////////////// float-double conversion ///////////////////////// |
1405 | 1453 | ||
1406 | int float2double(FVEC(x),DVEC(y)) { | 1454 | int float2double(FVEC(x),DVEC(y)) { |
@@ -1412,6 +1460,16 @@ int float2double(FVEC(x),DVEC(y)) { | |||
1412 | OK | 1460 | OK |
1413 | } | 1461 | } |
1414 | 1462 | ||
1463 | int float2int(KFVEC(x),IVEC(y)) { | ||
1464 | DEBUGMSG("float2int") | ||
1465 | int k; | ||
1466 | for(k=0;k<xn;k++) { | ||
1467 | yp[k]=xp[k]; | ||
1468 | } | ||
1469 | OK | ||
1470 | } | ||
1471 | |||
1472 | |||
1415 | int double2float(DVEC(x),FVEC(y)) { | 1473 | int double2float(DVEC(x),FVEC(y)) { |
1416 | DEBUGMSG("double2float") | 1474 | DEBUGMSG("double2float") |
1417 | int k; | 1475 | int k; |
@@ -1421,6 +1479,37 @@ int double2float(DVEC(x),FVEC(y)) { | |||
1421 | OK | 1479 | OK |
1422 | } | 1480 | } |
1423 | 1481 | ||
1482 | |||
1483 | int double2int(KDVEC(x),IVEC(y)) { | ||
1484 | DEBUGMSG("double2int") | ||
1485 | int k; | ||
1486 | for(k=0;k<xn;k++) { | ||
1487 | yp[k]=xp[k]; | ||
1488 | } | ||
1489 | OK | ||
1490 | } | ||
1491 | |||
1492 | |||
1493 | int int2float(KIVEC(x),FVEC(y)) { | ||
1494 | DEBUGMSG("int2float") | ||
1495 | int k; | ||
1496 | for(k=0;k<xn;k++) { | ||
1497 | yp[k]=xp[k]; | ||
1498 | } | ||
1499 | OK | ||
1500 | } | ||
1501 | |||
1502 | |||
1503 | int int2double(KIVEC(x),DVEC(y)) { | ||
1504 | DEBUGMSG("int2double") | ||
1505 | int k; | ||
1506 | for(k=0;k<xn;k++) { | ||
1507 | yp[k]=xp[k]; | ||
1508 | } | ||
1509 | OK | ||
1510 | } | ||
1511 | |||
1512 | |||
1424 | //////////////////// conjugate ///////////////////////// | 1513 | //////////////////// conjugate ///////////////////////// |
1425 | 1514 | ||
1426 | int conjugateQ(KQVEC(x),QVEC(t)) { | 1515 | int conjugateQ(KQVEC(x),QVEC(t)) { |
@@ -1465,8 +1554,30 @@ int stepD(DVEC(x),DVEC(y)) { | |||
1465 | OK | 1554 | OK |
1466 | } | 1555 | } |
1467 | 1556 | ||
1557 | |||
1468 | //////////////////// cond ///////////////////////// | 1558 | //////////////////// cond ///////////////////////// |
1469 | 1559 | ||
1560 | int compareF(KFVEC(x),KFVEC(y),IVEC(r)) { | ||
1561 | REQUIRES(xn==yn && xn==rn ,BAD_SIZE); | ||
1562 | DEBUGMSG("compareF") | ||
1563 | int k; | ||
1564 | for(k=0;k<xn;k++) { | ||
1565 | rp[k] = xp[k]<yp[k]?-1:(xp[k]>yp[k]?1:0); | ||
1566 | } | ||
1567 | OK | ||
1568 | } | ||
1569 | |||
1570 | int compareD(KDVEC(x),KDVEC(y),IVEC(r)) { | ||
1571 | REQUIRES(xn==yn && xn==rn ,BAD_SIZE); | ||
1572 | DEBUGMSG("compareD") | ||
1573 | int k; | ||
1574 | for(k=0;k<xn;k++) { | ||
1575 | rp[k] = xp[k]<yp[k]?-1:(xp[k]>yp[k]?1:0); | ||
1576 | } | ||
1577 | OK | ||
1578 | } | ||
1579 | |||
1580 | |||
1470 | int condF(FVEC(x),FVEC(y),FVEC(lt),FVEC(eq),FVEC(gt),FVEC(r)) { | 1581 | int condF(FVEC(x),FVEC(y),FVEC(lt),FVEC(eq),FVEC(gt),FVEC(r)) { |
1471 | REQUIRES(xn==yn && xn==ltn && xn==eqn && xn==gtn && xn==rn ,BAD_SIZE); | 1582 | REQUIRES(xn==yn && xn==ltn && xn==eqn && xn==gtn && xn==rn ,BAD_SIZE); |
1472 | DEBUGMSG("condF") | 1583 | DEBUGMSG("condF") |
@@ -1487,3 +1598,69 @@ int condD(DVEC(x),DVEC(y),DVEC(lt),DVEC(eq),DVEC(gt),DVEC(r)) { | |||
1487 | OK | 1598 | OK |
1488 | } | 1599 | } |
1489 | 1600 | ||
1601 | |||
1602 | int chooseF(KIVEC(cond),KFVEC(lt),KFVEC(eq),KFVEC(gt),FVEC(r)) { | ||
1603 | REQUIRES(condn==ltn && ltn==eqn && ltn==gtn && ltn==rn ,BAD_SIZE); | ||
1604 | DEBUGMSG("chooseF") | ||
1605 | int k; | ||
1606 | for(k=0;k<condn;k++) { | ||
1607 | rp[k] = condp[k]<0?ltp[k]:(condp[k]>0?gtp[k]:eqp[k]); | ||
1608 | } | ||
1609 | OK | ||
1610 | } | ||
1611 | |||
1612 | |||
1613 | int chooseD(KIVEC(cond),KDVEC(lt),KDVEC(eq),KDVEC(gt),DVEC(r)) { | ||
1614 | REQUIRES(condn==ltn && ltn==eqn && ltn==gtn && ltn==rn ,BAD_SIZE); | ||
1615 | DEBUGMSG("chooseD") | ||
1616 | int k; | ||
1617 | for(k=0;k<condn;k++) { | ||
1618 | rp[k] = condp[k]<0?ltp[k]:(condp[k]>0?gtp[k]:eqp[k]); | ||
1619 | } | ||
1620 | OK | ||
1621 | } | ||
1622 | |||
1623 | //////////////////////// extract ///////////////////////////////// | ||
1624 | |||
1625 | #define EXTRACT_IMP \ | ||
1626 | REQUIRES((tm == 0 && jn==rr && mc==rc) || (jn==rr && mr==rc) ,BAD_SIZE); \ | ||
1627 | DEBUGMSG("extractRD") \ | ||
1628 | int k,i,s; \ | ||
1629 | if (tm==0) { \ | ||
1630 | for (k=0;k<jn;k++) { \ | ||
1631 | s = jp[k]; \ | ||
1632 | for (i=0; i<mc; i++) { \ | ||
1633 | rp[rc*k+i] = mp[mc*s+i]; \ | ||
1634 | } \ | ||
1635 | } \ | ||
1636 | } else { \ | ||
1637 | for (k=0;k<jn;k++) { \ | ||
1638 | s = jp[k]; \ | ||
1639 | for (i=0; i<mr; i++) { \ | ||
1640 | rp[rc*k+i] = mp[mc*i+s]; \ | ||
1641 | } \ | ||
1642 | } \ | ||
1643 | } \ | ||
1644 | OK | ||
1645 | |||
1646 | |||
1647 | int extractRD(int tm, KIVEC(j), KDMAT(m), DMAT(r)) { | ||
1648 | EXTRACT_IMP | ||
1649 | } | ||
1650 | |||
1651 | int extractRF(int tm, KIVEC(j), KFMAT(m), FMAT(r)) { | ||
1652 | EXTRACT_IMP | ||
1653 | } | ||
1654 | |||
1655 | int extractRC(int tm, KIVEC(j), KCMAT(m), CMAT(r)) { | ||
1656 | EXTRACT_IMP | ||
1657 | } | ||
1658 | |||
1659 | int extractRQ(int tm, KIVEC(j), KQMAT(m), QMAT(r)) { | ||
1660 | EXTRACT_IMP | ||
1661 | } | ||
1662 | |||
1663 | int extractRI(int tm, KIVEC(j), KIMAT(m), IMAT(r)) { | ||
1664 | EXTRACT_IMP | ||
1665 | } | ||
1666 | |||