Python* API Reference for Intel® Data Analytics Acceleration Library 2019 Update 4

qr_dense_distr.py

Deprecation Notice: With the introduction of daal4py, a package that supersedes PyDAAL, Intel is deprecating PyDAAL and will discontinue support starting with Intel® DAAL 2021 and Intel® Distribution for Python 2021. Until then Intel will continue to provide compatible pyDAAL pip and conda packages for newer releases of Intel DAAL and make it available in open source. However, Intel will not add the new features of Intel DAAL to pyDAAL. Intel recommends developers switch to and use daal4py.

Note: To find daal4py examples, refer to daal4py documentation or browse github repository.

1 # file: qr_dense_distr.py
2 #===============================================================================
3 # Copyright 2014-2019 Intel Corporation.
4 #
5 # This software and the related documents are Intel copyrighted materials, and
6 # your use of them is governed by the express license under which they were
7 # provided to you (License). Unless the License provides otherwise, you may not
8 # use, modify, copy, publish, distribute, disclose or transmit this software or
9 # the related documents without Intel's prior written permission.
10 #
11 # This software and the related documents are provided as is, with no express
12 # or implied warranties, other than those that are expressly stated in the
13 # License.
14 #===============================================================================
15 
16 ## <a name="DAAL-EXAMPLE-PY-QR_DISTRIBUTED"></a>
17 ## \example qr_dense_distr.py
18 
19 import os
20 import sys
21 
22 from daal import step1Local, step2Master, step3Local
23 from daal.algorithms import qr
24 from daal.data_management import FileDataSource, DataSourceIface
25 
26 utils_folder = os.path.realpath(os.path.abspath(os.path.dirname(os.path.dirname(__file__))))
27 if utils_folder not in sys.path:
28  sys.path.insert(0, utils_folder)
29 from utils import printNumericTable
30 
31 DAAL_PREFIX = os.path.join('..', 'data')
32 
33 # Input data set parameters
34 nBlocks = 4
35 
36 datasetFileNames = [
37  os.path.join(DAAL_PREFIX, 'distributed', 'qr_1.csv'),
38  os.path.join(DAAL_PREFIX, 'distributed', 'qr_2.csv'),
39  os.path.join(DAAL_PREFIX, 'distributed', 'qr_3.csv'),
40  os.path.join(DAAL_PREFIX, 'distributed', 'qr_4.csv')
41 ]
42 
43 dataFromStep1ForStep2 = [0] * nBlocks
44 dataFromStep1ForStep3 = [0] * nBlocks
45 dataFromStep2ForStep3 = [0] * nBlocks
46 R = None
47 Qi = [0] * nBlocks
48 
49 
50 def computestep1Local(block):
51  global dataFromStep1ForStep2, dataFromStep1ForStep3
52 
53  # Initialize FileDataSource<CSVFeatureManager> to retrieve the input data from a .csv file
54  dataSource = FileDataSource(
55  datasetFileNames[block],
56  DataSourceIface.doAllocateNumericTable,
57  DataSourceIface.doDictionaryFromContext
58  )
59 
60  # Retrieve the input data
61  dataSource.loadDataBlock()
62 
63  # Create an algorithm to compute QR decomposition on the local node
64  algorithm = qr.Distributed(step1Local)
65 
66  algorithm.input.set(qr.data, dataSource.getNumericTable())
67 
68  # Compute QR decomposition and get OnlinePartialResult class from daal.algorithms.qr
69  pres = algorithm.compute()
70 
71  dataFromStep1ForStep2[block] = pres.get(qr.outputOfStep1ForStep2)
72  dataFromStep1ForStep3[block] = pres.get(qr.outputOfStep1ForStep3)
73 
74 
75 def computeOnMasterNode():
76  global R, dataFromStep2ForStep3
77 
78  # Create an algorithm to compute QR decomposition on the master node
79  algorithm = qr.Distributed(step2Master)
80 
81  for i in range(nBlocks):
82  algorithm.input.add(qr.inputOfStep2FromStep1, i, dataFromStep1ForStep2[i])
83 
84  # Compute QR decomposition and get DistributedPartialResult class from daal.algorithms.qr
85  pres = algorithm.compute()
86 
87  for i in range(nBlocks):
88  dataFromStep2ForStep3[i] = pres.getCollection(qr.outputOfStep2ForStep3, i)
89 
90  res = algorithm.finalizeCompute()
91  R = res.get(qr.matrixR)
92 
93 
94 def finalizeComputestep1Local(block):
95  global Qi
96 
97  # Create an algorithm to compute QR decomposition on the master node
98  algorithm = qr.Distributed(step3Local)
99 
100  algorithm.input.set(qr.inputOfStep3FromStep1, dataFromStep1ForStep3[block])
101  algorithm.input.set(qr.inputOfStep3FromStep2, dataFromStep2ForStep3[block])
102 
103  # Compute QR decomposition
104  algorithm.compute()
105 
106  res = algorithm.finalizeCompute()
107 
108  Qi[block] = res.get(qr.matrixQ)
109 
110 if __name__ == "__main__":
111 
112  for i in range(nBlocks):
113  computestep1Local(i)
114 
115  computeOnMasterNode()
116 
117  for i in range(nBlocks):
118  finalizeComputestep1Local(i)
119 
120  # Print the results
121  printNumericTable(Qi[0], "Part of orthogonal matrix Q from 1st node:", 10)
122  printNumericTable(R, "Triangular matrix R:")

For more complete information about compiler optimizations, see our Optimization Notice.