# # mimo.r # # Konstanten und Parameter setzen # Empfangsantennen (nicht ändern wegen Determinantenberechnung!) r <- 3 # Anzahl der Sendeantennen variabel bis tmax tmax <- 15 # Varianz der Kanalmatrix sigma <- 1 # Varianz des Rauschens N0 <- 1 # Varianz der Fehlermatrix fsigma <- sigma/10 # Leistungsbeschränkung P <- 10 # Anzahl der Durchläufe zur Durchschnittsermittlung n <- 10000 # Initialisierung der Kapazitäten/Schranken ev <- array(data = 0, dim = tmax) sv <- array(data = 0, dim = tmax) epu <- array(data = 0, dim = tmax) epo <- array(data = 0, dim = tmax) spu <- array(data = 0, dim = tmax) spo <- array(data = 0, dim = tmax) evk <- array(data = 0, dim = tmax) svk <- array(data = 0, dim = tmax) epuk <- array(data = 0, dim = tmax) epok <- array(data = 0, dim = tmax) spuk <- array(data = 0, dim = tmax) spok <- array(data = 0, dim = tmax) # Fehlertoleranz fehler <- 0.000001 Hmax <- array(data = 0, dim = r) Hi <- array(data = 0, dim = r) # Einheitsmatrix I <- diag(r) # Schleife: Anzahl der Durchläufe für Mittelung for (z in 1:n) { # Schleife: variable Anzahl von Sendeantennen for (t in 1:tmax) { H <- matrix(nrow = r, ncol = t) for (i in 1:r) { for (j in 1:t) { H[i,j] <- complex(real = rnorm(1, mean=0, sd=sigma), imag = rnorm(1, mean=0, sd=sigma)) } } HHstern <- H %*% Conj(t(H)) # Senderseitig vollständig bekannt lambda <- eigen(HHstern)$values my <- P/min(r,t) for (i in 1:min(r,t)) { my <- my + lambda[i]^-1/min(r,t) } summe <- 0 for (i in 1:min(r,t)) { if (my>lambda[i]^-1) summe <- summe + (my - lambda[i]^-1) } if (abs(P - summe) > fehler) { my <- P/(min(r,t)-1) for (i in 1:(min(r,t)-1)) { my <- my + lambda[i]^-1/(min(r,t)-1) } summe <- 0 for (i in 1:(min(r,t)-1)) { if (my>lambda[i]^-1) summe <- summe + (my - lambda[i]^-1) } if (abs(P - summe) > fehler) { # Dieser Fall tritt so selten ein, dass er nicht berücksichtigt wird. cat("Fehler! Bitte starten Sie die Simulation neu. \n") } } summe <- 0 for (i in 1:min(r,t)) { if (my>0) if (log(my*lambda[i], base=exp(1))>0) summe <- summe + log(my*lambda[i], base=exp(1)) } sv[t] <- sv[t] + summe # Empfängerseitig vollständig bekannt evmat <- I + P/t * HHstern evdet <- evmat[1,1]*evmat[2,2]*evmat[3,3] + evmat[2,1]*evmat[3,2]*evmat[1,3] + evmat[1,2]*evmat[2,3]*evmat[3,1] - evmat[1,3]*evmat[2,2]*evmat[3,1] - evmat[2,3]*evmat[3,2]*evmat[1,1] - evmat[1,2]*evmat[2,1]*evmat[3,3] # Imaginärteil = 0. Rundungsfehler müssen ausgeschlossen werden. evdet <- Re(evdet) ev[t] <- ev[t] + log(evdet, base=exp(1)) # Empfängerseitig partiell bekannt (untere Schranke) epumat <- I + P/t/(P*sigma+N0) * HHstern epudet <- epumat[1,1]*epumat[2,2]*epumat[3,3] + epumat[2,1]*epumat[3,2]*epumat[1,3] + epumat[1,2]*epumat[2,3]*epumat[3,1] - epumat[1,3]*epumat[2,2]*epumat[3,1] - epumat[2,3]*epumat[3,2]*epumat[1,1] - epumat[1,2]*epumat[2,1]*epumat[3,3] # Imaginärteil = 0. Rundungsfehler müssen ausgeschlossen werden. epudet <- Re(epudet) epu[t] <- epu[t] + log(epudet, base=exp(1)) # Empfängerseitig partiell bekannt (obere Schranke) epomat <- I + 1/sigma * HHstern epodet <- epomat[1,1]*epomat[2,2]*epomat[3,3] + epomat[2,1]*epomat[3,2]*epomat[1,3] + epomat[1,2]*epomat[2,3]*epomat[3,1] - epomat[1,3]*epomat[2,2]*epomat[3,1] - epomat[2,3]*epomat[3,2]*epomat[1,1] - epomat[1,2]*epomat[2,1]*epomat[3,3] # Imaginärteil = 0. Rundungsfehler müssen ausgeschlossen werden. epodet <- Re(epodet) epo[t] <- epo[t] + r*log(P*sigma/t/N0) + log(epodet, base=exp(1)) # Senderseitig partiell bekannt (untere Schranke) spumat <- I + P/t/(P*fsigma+N0) * HHstern spudet <- spumat[1,1]*spumat[2,2]*spumat[3,3] + spumat[2,1]*spumat[3,2]*spumat[1,3] + spumat[1,2]*spumat[2,3]*spumat[3,1] - spumat[1,3]*spumat[2,2]*spumat[3,1] - spumat[2,3]*spumat[3,2]*spumat[1,1] - spumat[1,2]*spumat[2,1]*spumat[3,3] # Imaginärteil = 0. Rundungsfehler müssen ausgeschlossen werden. spudet <- Re(spudet) spu[t] <- spu[t] + log(spudet, base=exp(1)) # Senderseitig partiell bekannt (obere Schranke) for (i in 1:r) { Hmax[i] <- max(abs(H[i,1:t])^2) summeR <- 0 if (t>1) for (l in 1:(t-1)) { summeR <- summeR + abs(H[i,l])^2*sum(abs(H[i,(l+1):t])^2) } Hi[i] <- sqrt(summeR) spo[t] <- spo[t] + log(P*(4*Hi[i]+Hmax[i]+fsigma)/N0+1, base=exp(1)) } # Hier folgen analoge Berechnungen für den Fall, dass P/t konstant gehalten wird - nicht P. # Dazu wird in den obigen Rechnungen P ersetzt durch t*P, so dass das fest definierte P die Leistungsbeschränkung pro Antenne angibt. # Senderseitig vollständig bekannt lambda <- eigen(HHstern)$values my <- P*t/min(r,t) for (i in 1:min(r,t)) { my <- my + lambda[i]^-1/min(r,t) } summe <- 0 for (i in 1:min(r,t)) { if (my>lambda[i]^-1) summe <- summe + (my - lambda[i]^-1) } if (abs(P*t - summe) > fehler) { my <- P*t/(min(r,t)-1) for (i in 1:(min(r,t)-1)) { my <- my + lambda[i]^-1/(min(r,t)-1) } summe <- 0 for (i in 1:(min(r,t)-1)) { if (my>lambda[i]^-1) summe <- summe + (my - lambda[i]^-1) } if (abs(P*t - summe) > fehler) { # Dieser Fall tritt so selten ein, dass er nicht berücksichtigt wird. print("Fehler! Bitte starten Sie die Simulation neu.") } } summe <- 0 for (i in 1:min(r,t)) { if (my>0) if (log(my*lambda[i], base=exp(1))>0) summe <- summe + log(my*lambda[i], base=exp(1)) } svk[t] <- svk[t] + summe # Empfängerseitig vollständig bekannt evmat <- I + P * HHstern evdet <- evmat[1,1]*evmat[2,2]*evmat[3,3] + evmat[2,1]*evmat[3,2]*evmat[1,3] + evmat[1,2]*evmat[2,3]*evmat[3,1] - evmat[1,3]*evmat[2,2]*evmat[3,1] - evmat[2,3]*evmat[3,2]*evmat[1,1] - evmat[1,2]*evmat[2,1]*evmat[3,3] # Imaginärteil = 0. Rundungsfehler müssen ausgeschlossen werden. evdet <- Re(evdet) evk[t] <- evk[t] + log(evdet, base=exp(1)) # Empfängerseitig partiell bekannt (untere Schranke) epumat <- I + P/(P*t*sigma+N0) * HHstern epudet <- epumat[1,1]*epumat[2,2]*epumat[3,3] + epumat[2,1]*epumat[3,2]*epumat[1,3] + epumat[1,2]*epumat[2,3]*epumat[3,1] - epumat[1,3]*epumat[2,2]*epumat[3,1] - epumat[2,3]*epumat[3,2]*epumat[1,1] - epumat[1,2]*epumat[2,1]*epumat[3,3] # Imaginärteil = 0. Rundungsfehler müssen ausgeschlossen werden. epudet <- Re(epudet) epuk[t] <- epuk[t] + log(epudet, base=exp(1)) # Empfängerseitig partiell bekannt (obere Schranke) epomat <- I + 1/sigma * HHstern epodet <- epomat[1,1]*epomat[2,2]*epomat[3,3] + epomat[2,1]*epomat[3,2]*epomat[1,3] + epomat[1,2]*epomat[2,3]*epomat[3,1] - epomat[1,3]*epomat[2,2]*epomat[3,1] - epomat[2,3]*epomat[3,2]*epomat[1,1] - epomat[1,2]*epomat[2,1]*epomat[3,3] # Imaginärteil = 0. Rundungsfehler müssen ausgeschlossen werden. epodet <- Re(epodet) epok[t] <- epok[t] + r*log(P*sigma/N0) + log(epodet, base=exp(1)) # Senderseitig partiell bekannt (untere Schranke) spumat <- I + P/(P*t*fsigma+N0) * HHstern spudet <- spumat[1,1]*spumat[2,2]*spumat[3,3] + spumat[2,1]*spumat[3,2]*spumat[1,3] + spumat[1,2]*spumat[2,3]*spumat[3,1] - spumat[1,3]*spumat[2,2]*spumat[3,1] - spumat[2,3]*spumat[3,2]*spumat[1,1] - spumat[1,2]*spumat[2,1]*spumat[3,3] # Imaginärteil = 0. Rundungsfehler müssen ausgeschlossen werden. spudet <- Re(spudet) spuk[t] <- spuk[t] + log(spudet, base=exp(1)) # Senderseitig partiell bekannt (obere Schranke) # Die Arrays Hmax und Hi wurden schon oben berechnet for (i in 1:r) { spok[t] <- spok[t] + log(P*t*(4*Hi[i]+Hmax[i]+fsigma)/N0+1, base=exp(1)) } } } sv <- sv/n ev <- ev/n epu <- epu/n epo <- epo/n spu <- spu/n spo <- spo/n svk <- svk/n evk <- evk/n epuk <- epuk/n epok <- epok/n spuk <- spuk/n spok <- spok/n x11(width = 7, height = 7, pointsize = 12) plot(1:15, 2*1:15, type="n", main="Verlauf der Kapazität bei wachsendem t und konst. P", xlab="t", ylab="C (Bits/Zeiteinheit)") lines(sv, col="blue", lty="solid") lines(ev, col="green", lty="dashed") lines(epu, col="black", lty="dotted") lines(epo, col="black", lty="dotted") lines(spu, col="red", lty="dotdash") lines(spo, col="red", lty="dotdash") x11(width = 7, height = 7, pointsize = 12) plot(1:15, 2*1:15, type="n", main="Verlauf der Kapazität bei wachsendem t und konst. P/t", xlab="t", ylab="C (Bits/Zeiteinheit)") lines(svk, col="blue", lty="solid") lines(evk, col="green", lty="dashed") lines(epuk, col="black", lty="dotted") lines(epok, col="black", lty="dotted") lines(spuk, col="red", lty="dotdash") lines(spok, col="red", lty="dotdash") x11(width = 7, height = 7, pointsize = 12) plot(1, 1, axes=FALSE, type="n", xlab="", ylab="", main="Legende") lines(c(0,1,2), c(1.4,1.4,1.4), col="blue", lty="solid") lines(c(0,1,2), c(1.2,1.2,1.2), col="green", lty="dashed") lines(c(0,1,2), c(1,1,1), col="black", lty="dotted") lines(c(0,1,2), c(0.8,0.8,0.8), col="red", lty="dotdash") text(1, 1.3, "senderseitig vollständig bekannt") text(1, 1.1, "empfängerseitig vollständig bekannt") text(1, 0.9, "empfängerseitig partiell bekannt (obere und untere Schranke)") text(1, 0.7, "senderseitig partiell bekannt (obere und untere Schranke)")