summaryrefslogtreecommitdiff
path: root/2.3-1/src/c/linearAlgebra/qr/dqra.c
diff options
context:
space:
mode:
Diffstat (limited to '2.3-1/src/c/linearAlgebra/qr/dqra.c')
-rw-r--r--2.3-1/src/c/linearAlgebra/qr/dqra.c298
1 files changed, 0 insertions, 298 deletions
diff --git a/2.3-1/src/c/linearAlgebra/qr/dqra.c b/2.3-1/src/c/linearAlgebra/qr/dqra.c
deleted file mode 100644
index bae4bc27..00000000
--- a/2.3-1/src/c/linearAlgebra/qr/dqra.c
+++ /dev/null
@@ -1,298 +0,0 @@
-/* Copyright (C) 2017 - IIT Bombay - FOSSEE
-
- This file must be used under the terms of the CeCILL.
- This source file is licensed as described in the file COPYING, which
- you should have received as part of this distribution. The terms
- are also available at
- http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
- Author: Sandeep Gupta
- Organization: FOSSEE, IIT Bombay
- Email: toolbox@scilab.in
- */
-
-/* This C code is used to generate function for QR decomposition */
-
-#include <stdio.h>
-#include <stdlib.h>
-#include <math.h>
-#include "qr.h"
-#include "lapack.h"
-#include "string.h"
-#include "matrixTranspose.h"
-
-/*For reference check Scilab source code & lapack library websites
-Names of variable are almost same for convience.
-*/
-/*
-int min(int M,int N){
- if(M > N) return N;
- return M;
-}
-
-int max(int M,int N){
- if(M > N) return N;
- return N;
-}*/
-
-/* External Function used of lapack library */
-extern double dgeqrf_(int *,int *,double *,int *,double *,double *,int *,int *);
-extern double dlacpy_(char *,int *,int *,double *,int *,double *,int *);
-extern double dgeqpf_(int *,int *,double *,int *,int *,double *,double *,int *);
-extern double dorgqr_(int *,int *,int *,double *,int *,double *,double *,int *,int *);
-extern void dlaset_(char *,int *,int *,double *,double *,double *,int *);
-
-/* function for finding qr */
-double dqra(int ninp,int nout,double *inp1,int M,int N,double tol,double *out1,double *out2,double *out3){
- int i,j;
- char choice;
- double alpha=0.0,beta=0.0;
- int minMN = min(M,N);
-
- double *A;
- A = (double *)malloc(M*N*sizeof(double));
- memcpy(A,inp1,M*N*sizeof(double));
-
- if(M <= 0 || N <= 0){
- out1 = NULL;
- out2 = NULL;
- return 0;
- }
- /* doldqr */
- if(nout == 4){ /* [Q,R,rk,E]=qr(X [,tol]) */
- if(ninp == 1){
- tol = -1;
- }
-
- int INFO,rk;
-
- int *JPVT;
- JPVT = ( int *)malloc(N*sizeof(int));
-
- double *TAU;
- TAU = (double *)malloc(minMN*sizeof(double));
-
- int LWORK = 3*N;
-
- double *WORK;
- WORK = (double *)malloc(LWORK*sizeof(double));
-
- for(i=1;i<=M;i++){
- JPVT[i-1]=0.0;
- }
- dgeqpf_(&M,&N,A,&M,JPVT,TAU,WORK,&INFO);
-
- choice = 'U';
- dlacpy_(&choice,&M,&N,A,&M,out2,&M);
-
- if(M > N){
- for(j=1;j<=N;j++){
- for(i=j+1;i<=M;i++){
- out2[i-1+(j-1)*M] = 0.0;
- }
- }
- }
- else{
- for(j=1;j<=M-1;j++){
- for(i=j+1;i<=M;i++){
- out2[i-1+(j-1)*M] = 0.0;
- }
- }
- }
-
- if(M > N){
- choice = 'F';
- dlacpy_(&choice,&M,&N,A,&M,out1,&M);
- for(j=N+1;j<=M;j++){
- for(i=1;i<=M;i++){
- out1[i-1+(j-1)*M] = 0.0;
- }
- }
- }
- else{
- choice = 'F';
- dlacpy_(&choice,&M,&M,A,&M,out1,&M);
- }
-
- dorgqr_(&M,&M,&minMN,out1,&M,TAU,WORK,&LWORK,&INFO);
-
- choice = 'F';
- dlaset_(&choice,&N,&N,&alpha,&beta,out3,&N);
-
- for(j=1;j<=N;j++){
- i = JPVT[j-1];
- out3[i-1+(j-1)*N] = 1.0;
- }
-
- double tt = abs(out2[0]);
-
- if(tol == -1){
- tol = (double)max(M,N)*pow(2,-52)*tt;
- }
- double ch;
- //printf("%lf ",tol);
- for(j=1;j<=minMN;j++){
- //printf("%e ",out2[(j-1)+(j-1)*M]);
- if(out2[(j-1)+(j-1)*M] < 0){
- ch = -out2[(j-1)+(j-1)*M];
- }
- else{
- ch = out2[(j-1)+(j-1)*M];
- }
- if(ch > tol){
- rk = j;
- }
- else{
- break;
- }
- }
- //printf("\n");
- return rk;
- }
- else if(ninp == 1){
- /* (intdgeqpf3)
- [Q,R]=qr(A)
- [Q,R,E]=qr(A)
- */
- int LDA = M;
-
- double *TAU;
- TAU = (double *)malloc(min(M,N)*sizeof(double));
-
- int LWORK;
- if(nout <= 2)
- LWORK = N;
- else
- LWORK = 3*N;
-
- double *WORK;
- WORK = (double *)malloc((LWORK+1)*sizeof(double));
-
- int INFO;
-
- int *JPVT;
- JPVT = (int *)malloc(N*sizeof(int));
-
- if(nout <= 2){
- dgeqrf_(&M,&N,A,&M,TAU,WORK,&LWORK,&INFO);
- }
- else{
- for(i=0;i<N;i++){
- *(JPVT+i)=0.0;
- }
- dgeqpf_(&M,&N,A,&M,JPVT,TAU,WORK,&INFO);
- }
-
- /* Copying code from A to R */
- choice = 'U';
- dlacpy_(&choice,&M,&N,A,&M,out2,&M);
- /*for(i=0;i<M;i++){
- for(j=0;j<N-i;j++){
- out2[i+j*M] = A[i+j*M];
- }
- }*/
-
- if(M > N){
- for(j=1;j<=N;j++){
- for(i=j+1;i<=M;i++){
- out2[i-1+(j-1)*M] = 0.0;
- }
- }
- }
- else{
- for(j=1;j<=M-1;j++){
- for(i=j+1;i<=M;i++){
- out2[i-1+(j-1)*M] = 0.0;
- }
- }
- }
- /* lQ - out1 */
- if(M > N){
- choice = 'F';
- dlacpy_(&choice,&M,&N,A,&M,out1,&M);
- /*for(i=0;i<M*N;i++){
- out1[i]=A[i];
- }*/
- for(j=N+1;j<=M;j++){
- for(i=1;i<=M;i++){
- out1[i+(j-1)*M] = 0.0;
- }
- }
- }
- else{
- choice = 'F';
- dlacpy_(&choice,&M,&M,A,&M,out1,&M);
- /*for(i=0;i<M;i++){
- for(j=0;j<M;j++){
- out1
- }
- }*/
- }
- dorgqr_(&M,&M,&minMN,out1,&M,TAU,WORK,&LWORK,&INFO);
-
- if(nout > 2){
- choice = 'F';
- dlaset_(&choice,&N,&N,&alpha,&beta,out3,&N);
- for(j=1;j<=N;j++){
- i = *(JPVT+j-1);
- //printf("%d ",i-1+(j-1)*N);
- *(out3+i-1+(j-1)*N) = 1.0;
- }
- //printf("\n");
- }
- }
- else{/*
- [[Q,R]=qr(A,'e')
- [Q,R,E]=qr(A,'e') ] */
- int *JPVT;
- JPVT = (int *)malloc(N*sizeof(int));
-
- double *TAU;
- TAU = (double *)malloc(min(M,N)*sizeof(double));
- int LWORK;
- if(nout <= 2){
- LWORK = N;
- }
- else{
- LWORK = 3*N;
- }
-
- double *WORK;
- WORK = (double *)malloc(LWORK*sizeof(double));
-
- int INFO;
- if(nout <= 2)
- dgeqrf_(&M,&N,A,&M,TAU,WORK,&LWORK,&INFO);
- else{
- for(i=1;i<=N;i++){
- JPVT[i-1]=0.0;
- }
- dgeqpf_(&M,&N,A,&M,JPVT,TAU,WORK,&INFO);
- }
- choice = 'U';
-
- dlacpy_(&choice,&minMN,&N,A,&M,out2,&minMN);
-
- if(N >= 2){
- for(j=1;j<=N-1;j++){
- if(j+1 <= minMN){
- for(i=j+1;i<=minMN;i++){
- out2[i-1+(j-1)*minMN] = 0.0;
- }
- }
- }
- }
- choice = 'F';
- dlacpy_(&choice,&M,&minMN,A,&M,out1,&M);
- dorgqr_(&M,&minMN,&minMN,out1,&M,TAU,WORK,&LWORK,&INFO);
-
- if(nout > 2){
- choice = 'F';
- dlaset_(&choice,&N,&N,&alpha,&beta,out3,&N);
- for(j=1;j<=N;j++){
- i = JPVT[j-1];
- out3[i-1+(j-1)*N] = 1.0;
- }
- }
- }
- return 0;
-}